diff --git a/src/geometry/structured/BinaryRawFile.cpp b/src/geometry/structured/BinaryRawFile.cpp new file mode 100644 index 0000000000000000000000000000000000000000..5b4c81256db3041cc6617fd245720f421da70915 --- /dev/null +++ b/src/geometry/structured/BinaryRawFile.cpp @@ -0,0 +1,71 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file BinaryRawFileReader.cpp +//! \ingroup geometry +//! \author Christian Godenschwager <christian.godenschwager@fau.de> +// +//====================================================================================================================== + +#include "BinaryRawFile.h" + + +namespace walberla { +namespace geometry { + +BinaryRawFile::BinaryRawFile( const std::string & filename, const Vector3< uint_t > & size, const std::string & datatype ) + : size_( size ) +{ + init( filename, datatype ); +} + +BinaryRawFile::BinaryRawFile( const Config::BlockHandle & configBlock ) +{ + size_ = configBlock.getParameter< Vector3< uint_t > >( "size" ); + const std::string filename = configBlock.getParameter< std::string >( "filename" ); + const std::string datatype = configBlock.getParameter< std::string >( "datatype" ); + + init( filename, datatype ); +} + +void BinaryRawFile::init( const std::string & filename, const std::string & datatype ) +{ + if (datatype == "uint8") + init<uint8_t>( filename ); + else if (datatype == "int8") + init<int8_t>( filename ); + else if (datatype == "uint16") + init<uint16_t>( filename ); + else if (datatype == "int16") + init<int16_t>( filename ); + else if (datatype == "uint32") + init<uint32_t>( filename ); + else if (datatype == "int32") + init<int32_t>( filename ); + else if (datatype == "uint64") + init<uint64_t>( filename ); + else if (datatype == "int64") + init<int64_t>( filename ); + else if (datatype == "float") + init<float>( filename ); + else if (datatype == "double") + init<double>( filename ); + else + WALBERLA_ABORT( "Unknown datatype!" ); +} + + +} // namespace geometry +} // namespace walberla \ No newline at end of file diff --git a/src/geometry/structured/BinaryRawFile.h b/src/geometry/structured/BinaryRawFile.h new file mode 100644 index 0000000000000000000000000000000000000000..8a9e69346fbc801860146c493ed5e7451d35d809 --- /dev/null +++ b/src/geometry/structured/BinaryRawFile.h @@ -0,0 +1,170 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file BinaryRawFileReader.h +//! \ingroup geometry +//! \author Christian Godenschwager <christian.godenschwager@fau.de> +// +//====================================================================================================================== + +#pragma once + +#include "core/DataTypes.h" +#include "core/mpi/Broadcast.h" +#include "core/mpi/MPIManager.h" +#include "core/mpi/BufferDataTypeExtensions.h" +#include "core/math/AABB.h" +#include "core/math/Vector3.h" +#include "core/config/Config.h" + +#include <fstream> +#include <string> +#include <vector> +#include <iterator> + + +namespace walberla { +namespace geometry { + +class BinaryRawFile +{ +public: + BinaryRawFile( const std::string & filename, const Vector3< uint_t > & size, const std::string & datatype ); + + template< typename T > + static BinaryRawFile loadFile( const std::string & filename, const Vector3< uint_t > & size ); + + BinaryRawFile( const Config::BlockHandle & configBlock ); + + inline bool get( const Vector3< uint_t > & pos ) const; + inline bool get( const uint_t x, const uint_t y, const uint_t z ) const; + + const Vector3< uint_t> & size() const { return size_; } + +private: + template< typename T > + BinaryRawFile( const std::string & filename, const Vector3< uint_t > & size, const T ); + + void init( const std::string & filename, const std::string & datatype ); + + template< typename T > + void init( const std::string & filename ); + + Vector3< uint_t > size_; + std::vector< bool > data_; +}; + +class BinaryRawFileInterpolator +{ +public: + enum Interpolator { NEAREST_NEIGHBOR }; + + BinaryRawFileInterpolator( const AABB & aabb, const BinaryRawFile & binaryRawFile, Interpolator interpolator ) + : aabb_(aabb), interpolator_(interpolator), binaryRawFile_( binaryRawFile ) {} + + inline bool get( const Vector3< real_t > & pos ) const; + inline bool get( const real_t x, const real_t y, const real_t z ) const; + + const AABB & aabb() const { return aabb_; } + + Interpolator interpolator() const { return interpolator_; } + +private: + inline bool getNearestNeighbor( const real_t x, const real_t y, const real_t z ) const; + + AABB aabb_; + Interpolator interpolator_; + const BinaryRawFile & binaryRawFile_; +}; + + +template<typename T> +BinaryRawFile::BinaryRawFile( const std::string & filename, const Vector3<uint_t> & size, const T ) + : size_(size) +{ + init<T>( filename ); +} + + +template<typename T> +void BinaryRawFile::init( const std::string & filename ) +{ + const uint_t numElements = size_[0] * size_[1] * size_[2]; + data_.reserve( numElements ); + WALBERLA_ROOT_SECTION() + { + std::ifstream ifs( filename, std::ifstream::in | std::ifstream::binary ); + std::transform( std::istream_iterator<T>( ifs ), std::istream_iterator<T>(), + std::back_inserter( data_ ), []( const T v ) { return v > T( 0 ); } ); + WALBERLA_CHECK_EQUAL( numElements, data_.size(), "Error reading file \"" << filename << "\"!" ); + } + mpi::broadcastObject( data_ ); +} + +template< typename T > +static BinaryRawFile loadFile( const std::string & filename, const Vector3< uint_t > & size ) +{ + return BinaryRawFileReader( filename, size, T() ); +} + +bool BinaryRawFile::get( const Vector3< uint_t > & pos ) const +{ + return get( pos[0], pos[1], pos[2] ); +} + +inline bool BinaryRawFile::get( const uint_t x, const uint_t y, const uint_t z ) const +{ + WALBERLA_ASSERT_LESS( x, size_[0] ); + WALBERLA_ASSERT_LESS( y, size_[1] ); + WALBERLA_ASSERT_LESS( z, size_[2] ); + + const uint_t i = z * size_[0] * size_[1] + y * size_[0] + x; + + WALBERLA_ASSERT_LESS( i, data_.size() ); + + return data_[i]; +} + + +bool BinaryRawFileInterpolator::get( const Vector3< real_t > & pos ) const +{ + return get( pos[0], pos[1], pos[2] ); +} + + +bool BinaryRawFileInterpolator::get( const real_t x, const real_t y, const real_t z ) const +{ + WALBERLA_ASSERT( aabb_.contains( x, y, z ) ); + + switch (interpolator_) + { + case NEAREST_NEIGHBOR: + return getNearestNeighbor( x, y, z ); + default: + WALBERLA_ABORT( "Unknown Interpolator!" ); + } +} + +bool BinaryRawFileInterpolator::getNearestNeighbor( const real_t x, const real_t y, const real_t z ) const +{ + uint_t xInt = uint_c( (x - aabb_.xMin()) / aabb_.xSize() * real_t( binaryRawFile_.size()[0] ) + real_t( 0.5 ) ); + uint_t yInt = uint_c( (y - aabb_.yMin()) / aabb_.ySize() * real_t( binaryRawFile_.size()[1] ) + real_t( 0.5 ) ); + uint_t zInt = uint_c( (z - aabb_.zMin()) / aabb_.zSize() * real_t( binaryRawFile_.size()[2] ) + real_t( 0.5 ) ); + + return binaryRawFile_.get( xInt, yInt, zInt ); +} + +} // namespace geometry +} // namespace walberla \ No newline at end of file diff --git a/tests/geometry/BinaryRawFileTest.cpp b/tests/geometry/BinaryRawFileTest.cpp new file mode 100644 index 0000000000000000000000000000000000000000..cc219ea1b162c98abaa93f246f6a3c841f276c5f --- /dev/null +++ b/tests/geometry/BinaryRawFileTest.cpp @@ -0,0 +1,155 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file BinaryRawFileTest.cpp +//! \ingroup geometry +//! \author Christian Godenschwager <christian.godenschwager@fau.de> +// +//====================================================================================================================== + +#include "blockforest/Initialization.h" + +#include "core/DataTypes.h" +#include "core/Abort.h" +#include "core/logging/Logging.h" +#include "core/mpi/Environment.h" + +#include "field/AddToStorage.h" +#include "field/vtk/VTKWriter.h" + +#include "geometry/structured/BinaryRawFile.h" + +#include "vtk/Initialization.h" + +#include <string> +#include <vector> + +namespace walberla { +namespace geometry { + +void test( const std::string & filename, const Vector3<uint_t> & size, const std::string & datatype, const bool vtkOut ) +{ + WALBERLA_LOG_INFO( "Testing unscaled" ); + BinaryRawFile brf( filename, size, datatype ); + + auto blocks = blockforest::createUniformBlockGrid( uint_t(1), uint_t( 1 ), uint_t( 1 ), + size[0], size[1], size[2], + real_t( 1 ), + uint_t( 1 ), uint_t( 1 ), uint_t( 1 ) ); + + typedef GhostLayerField< uint8_t, uint_t( 1 ) > ScalarField; + + BlockDataID scalarFieldID = field::addToStorage<ScalarField>( blocks, "BinaryRawFile" ); + + for ( auto & block : *blocks) + { + auto field = block.getData<ScalarField>( scalarFieldID ); + + CellInterval ci( 0, 0, 0, cell_idx_c( size[0] ) - 1, cell_idx_c( size[1] ) - 1, cell_idx_c( size[2] ) - 1 ); + + for (const Cell & c : ci) + { + field->get( c ) = brf.get( uint_c( c[0] ), uint_c( c[1] ), uint_c( c[2] ) ); + } + } + + if (vtkOut) + { + WALBERLA_LOG_INFO( "Writing unscaled" ); + auto vtkOutput = vtk::createVTKOutput_BlockData( blocks, "BinaryRawFile" ); + vtkOutput->addCellDataWriter( make_shared< field::VTKWriter< ScalarField, uint8_t > >( scalarFieldID, "BinaryRawFile" ) ); + writeFiles( vtkOutput, true )(); + } +} + +void testScaled( const std::string & filename, const Vector3<uint_t> & size, const std::string & datatype, const bool vtkOut ) +{ + WALBERLA_LOG_INFO( "Testing scaled" ); + BinaryRawFile brf( filename, size, datatype ); + + Vector3<uint_t> scaledSize( std::max( uint_t( 1 ), size[0] / uint_t( 2 ) ), + std::max( uint_t( 1 ), size[1] / uint_t( 3 ) ), + std::max( uint_t( 1 ), size[2] / uint_t( 5 ) ) ); + + auto blocks = blockforest::createUniformBlockGrid( uint_t( 1 ), uint_t( 1 ), uint_t( 1 ), + scaledSize[0], scaledSize[1], scaledSize[2], + real_t( 1 ), + uint_t( 1 ), uint_t( 1 ), uint_t( 1 ) ); + + BinaryRawFileInterpolator brfi( blocks->getDomain(), brf, BinaryRawFileInterpolator::NEAREST_NEIGHBOR ); + + typedef GhostLayerField< uint8_t, uint_t( 1 ) > ScalarField; + + BlockDataID scalarFieldID = field::addToStorage<ScalarField>( blocks, "BinaryRawFile" ); + + for (auto & block : *blocks) + { + auto field = block.getData<ScalarField>( scalarFieldID ); + + CellInterval ci( 0, 0, 0, cell_idx_c( scaledSize[0] ) - 1, cell_idx_c( scaledSize[1] ) - 1, cell_idx_c( scaledSize[2] ) - 1 ); + + for (const Cell & c : ci) + { + auto pos = blocks->getBlockLocalCellCenter( block, c ); + field->get( c ) = brfi.get( pos ); + } + } + + if (vtkOut) + { + WALBERLA_LOG_INFO( "Writing scaled" ); + auto vtkOutput = vtk::createVTKOutput_BlockData( blocks, "BinaryRawFileScaled" ); + vtkOutput->addCellDataWriter( make_shared< field::VTKWriter< ScalarField, uint8_t > >( scalarFieldID, "BinaryRawFile" ) ); + writeFiles( vtkOutput, true )(); + } +} + +} // namespace geometry +} // namespace walberla + +int main( int argc, char * argv[] ) +{ + walberla::mpi::Environment env( argc, argv ); + + std::vector< std::string > args( argv, argv + argc ); + + auto it = std::find( args.begin(), args.end(), std::string( "--vtk" ) ); + const bool vtk = it != args.end(); + if (it != args.end()) + args.erase( it ); + + if (args.size() != 6u) + WALBERLA_ABORT_NO_DEBUG_INFO( "USAGE: " << args.front() << " [--vtk] FILENAME XSIZE YSIZE ZSIZE DATATYPE" ); + + const std::string filename = args[1]; + + walberla::Vector3< walberla::uint_t > size; + try + { + size.set( std::stoull( args[2] ), + std::stoull( args[3] ), + std::stoull( args[4] ) ); + } + catch (const std::invalid_argument &) + { + WALBERLA_ABORT_NO_DEBUG_INFO( "USAGE: " << args.front() << " FILENAME XSIZE YSIZE ZSIZE DATATYPE" ); + } + + const std::string datatype = args[5]; + + walberla::geometry::test( filename, size, datatype, vtk ); + walberla::geometry::testScaled( filename, size, datatype, vtk ); +} + diff --git a/tests/geometry/CMakeLists.txt b/tests/geometry/CMakeLists.txt index 789e44395fb2012327da3299c2d1698dd35b65d4..fa297777fb0c149ca2f3ebdc7f8eea781c69f302 100644 --- a/tests/geometry/CMakeLists.txt +++ b/tests/geometry/CMakeLists.txt @@ -28,3 +28,5 @@ waLBerla_execute_test( NAME ScalarFieldFromGrayScaleImageTest ) waLBerla_compile_test( FILES ImageTest.cpp ) waLBerla_execute_test( NAME ImageTest ) + +waLBerla_compile_test( FILES BinaryRawFileTest.cpp )