songtianlun/openfoam
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

ENH: add support for additional filter/mapping (#2609)

Browse Source 
- comprises a few different elements:

FilterField (currently packaged in PatchFunction1Types namespace)
~~~~~~~~~~~

  The FilterField helper class provides a multi-sweep median filter
  for a Field of data associated with a geometric point cloud.

  The points can be freestanding or the faceCentres (or points)
  of a meshedSurface, for example.

  Using an initial specified search radius, the nearest point
  neighbours are gathered and addressing/weights are built for them.
  This currently uses an area-weighted, linear RBF interpolator
  with provision for quadratic RBF interpolator etc.

  After the weights and addressing are established,
  the evaluate() method can be called to apply a median filter
  to data fields, with a specified number of sweeps.

boundaryDataSurfaceReader
~~~~~~~~~~~~~~~~~~~~~~~~~

- a surfaceReader (similar to ensightSurfaceReader) when a general
  point data reader is needed.

MappedFile
~~~~~~~~~~
- has been extended to support alternative surface reading formats.
  This allows, for example, sampled ensight data to be reused for
  mapping.  Cavaet: multi-patch entries may still needs some work.

- additional multi-sweep median filtering of the input data.
  This can be used to remove higher spatial frequencies when
  sampling onto a coarse mesh.

smoothSurfaceData
~~~~~~~~~~~~~~~~~
- standalone application for testing of filter radii/sweeps
This commit is contained in:
Mark Olesen 2022-10-26 15:13:48 +02:00
parent 2984d1e3e7
commit cb4e026aed
13 changed files with 1790 additions and 8 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

3
applications/utilities/preProcessing/smoothSurfaceData/Make/files Normal file
View File

@ -0,0 +1,3 @@
smoothSurfaceData.C
EXE = $(FOAM_USER_APPBIN)/smoothSurfaceData

11
applications/utilities/preProcessing/smoothSurfaceData/Make/options Normal file
View File

@ -0,0 +1,11 @@
/* Can compile with -fopenmp etc */
EXE_INC = \
-I$(LIB_SRC)/fileFormats/lnInclude \
-I$(LIB_SRC)/surfMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lfileFormats \
-lsurfMesh \
-lmeshTools

180
applications/utilities/preProcessing/smoothSurfaceData/smoothSurfaceData.C Normal file
View File

@ -0,0 +1,180 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM 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.
OpenFOAM 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 OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Application
smoothSurfaceData
Description
Pre-processing, filtering of surface field data.
Currently this is just used to test filter settings
(for MappedFile) and only generates vtk output, which can be
easily loaded in ParaView.
\*---------------------------------------------------------------------------*/
#include "argList.H"
#include "Time.H"
#include "clockTime.H"
#include "primitiveFields.H"
#include "surfaceReader.H"
#include "surfaceWriter.H"
#include "foamVtkSurfaceWriter.H"
#include "MappedFileFilterField.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
using namespace Foam;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
argList::addNote
(
"Testing, pre-processing, filtering of surface field data"
);
argList::noCheckProcessorDirectories();
argList::addVerboseOption("Additional verbosity");
argList::addArgument("input", "The input surface file");
argList::addOption("radius", "m", "Specify filter radius [metres]");
argList::addOption("sweeps", "N", "Number of median filter stages");
argList::addOption
(
"field",
"name",
"Field <scalar> to process (default: T)"
);
argList::addOption
(
"read-format",
"type",
"Input format (default: ensight)"
);
argList args(argc, argv);
// #include "setRootCase.H"
const int optVerbose = args.verbose();
const word readFileType
(
args.getOrDefault<word>("read-format", "ensight")
);
// Constant radius searching
label filterSweeps_(1);
scalar filterRadius_(0);
args.readIfPresent("sweeps", filterSweeps_);
args.readIfPresent("radius", filterRadius_);
Info<< nl
<< "Filter: radius=" << filterRadius_
<< " sweeps=" << filterSweeps_ << endl;
// Simple sanity check
if ((filterSweeps_ < 1) || (filterRadius_ <= VSMALL))
{
Info<< nl << "Nothing to do. Exiting..." << nl << endl;
return 0;
}
word fieldName("T");
args.readIfPresent("field", fieldName);
const fileName importName = args.get<fileName>(1);
auto readerPtr_ = surfaceReader::New(readFileType, importName);
auto& reader = readerPtr_();
const label fieldIndex = reader.fieldNames(0).find(fieldName);
if (fieldIndex == -1)
{
FatalErrorInFunction
<< "Unable to find field name: " << fieldName
<< " in list of available fields: " << reader.fieldNames(0)
<< exit(FatalError);
}
clockTime timing;
const meshedSurface& geom = reader.geometry(0);
Info<< nl << "Read " << geom.nFaces() << " faces and "
<< geom.nPoints() << " points in "
<< timing.timeIncrement() << "s" << endl;
PatchFunction1Types::FilterField fieldFilter;
PatchFunction1Types::FilterField::debug = optVerbose;
fieldFilter.reset(geom, filterRadius_);
Info<< nl << "Built weights/addressing "
<< timing.timeIncrement() << "s" << endl;
Info<< nl << "Processing " << reader.times().size() << " times" << nl;
const instantList times(reader.times());
forAll(times, timeIndex)
{
tmp<scalarField> tfield
(
reader.field(timeIndex, fieldIndex, pTraits<scalar>::zero)
);
tfield = fieldFilter.evaluate(tfield, filterSweeps_);
vtk::surfaceWriter writer
(
geom.points(),
geom,
word::printf("filtered_%06d", timeIndex)
);
writer.piece(geom.points(), geom);
writer.writeGeometry();
writer.beginCellData();
writer.writeCellData(fieldName, tfield());
}
Info<< nl << "Smoothing/writing "
<< timing.timeIncrement() << "s" << endl;
Info<< nl << "End\n" << endl;
return 0;
}
// ************************************************************************* //

2
src/meshTools/Make/files
View File

@ -327,6 +327,8 @@ PatchFunction1/makePatchFunction1s.C
PatchFunction1/coordinateScaling/coordinateScalings.C
PatchFunction1/CodedField/makeCodedFields.C
PatchFunction1/MappedFile/MappedFileFilterField.C
meshStructure/meshStructure.C
coupling/externalFileCoupler.C

229
src/meshTools/PatchFunction1/MappedFile/MappedFile.C
View File

@ -47,6 +47,12 @@ Foam::PatchFunction1Types::MappedFile<Type>::MappedFile
fieldTableName_(dict.getOrDefault<word>("fieldTable", entryName)),
pointsName_(dict.getOrDefault<word>("points", "points")),
mapMethod_(),
filterRadius_(dict.getOrDefault<scalar>("filterRadius", 0)),
filterSweeps_(dict.getOrDefault<label>("filterSweeps", 0)),
filterFieldPtr_(nullptr),
readerFormat_(),
readerFile_(),
readerPtr_(nullptr),
mapperPtr_(nullptr),
sampleTimes_(),
sampleIndex_(-1, -1),
@ -54,6 +60,40 @@ Foam::PatchFunction1Types::MappedFile<Type>::MappedFile
sampleValues_(),
offset_(Function1<Type>::NewIfPresent("offset", dict))
{
// Simple sanity check
if ((filterSweeps_ < 1) || (filterRadius_ <= VSMALL))
{
filterRadius_ = 0;
filterSweeps_ = 0;
}
if (dict.readIfPresent("sampleFormat", readerFormat_))
{
dict.readEntry("sampleFile", readerFile_);
fileName fName(readerFile_);
fName.expand();
readerPtr_ = surfaceReader::New(readerFormat_, fName);
}
if (debug)
{
Info<< "mappedFile:" << nl;
if (readerFormat_.empty())
{
Info<< " boundary format" << nl;
}
else
{
Info<< " format:" << readerFormat_
<< " file:" << readerFile_ << nl;
}
Info<< " filter radius=" << filterRadius_
<< " sweeps=" << filterSweeps_ << endl;
}
if
(
dict.readIfPresent("mapMethod", mapMethod_)
@ -88,6 +128,12 @@ Foam::PatchFunction1Types::MappedFile<Type>::MappedFile
fieldTableName_(fieldTableName),
pointsName_(dict.getOrDefault<word>("points", "points")),
mapMethod_(),
filterRadius_(dict.getOrDefault<scalar>("filterRadius", 0)),
filterSweeps_(dict.getOrDefault<label>("filterSweeps", 0)),
filterFieldPtr_(nullptr),
readerFormat_(),
readerFile_(),
readerPtr_(nullptr),
mapperPtr_(nullptr),
sampleTimes_(),
sampleIndex_(-1, -1),
@ -95,6 +141,40 @@ Foam::PatchFunction1Types::MappedFile<Type>::MappedFile
sampleValues_(),
offset_(Function1<Type>::NewIfPresent("offset", dict))
{
// Simple sanity check
if ((filterSweeps_ < 1) || (filterRadius_ <= VSMALL))
{
filterRadius_ = 0;
filterSweeps_ = 0;
}
if (dict.readIfPresent("sampleFormat", readerFormat_))
{
dict.readEntry("sampleFile", readerFile_);
fileName fName(readerFile_);
fName.expand();
readerPtr_ = surfaceReader::New(readerFormat_, fName);
}
if (debug)
{
Info<< "mappedFile:" << nl;
if (readerFormat_.empty())
{
Info<< " boundary format" << nl;
}
else
{
Info<< " format:" << readerFormat_
<< " file:" << readerFile_ << nl;
}
Info<< " filter radius=" << filterRadius_
<< " sweeps=" << filterSweeps_ << endl;
}
if
(
dict.readIfPresent("mapMethod", mapMethod_)
@ -136,13 +216,27 @@ Foam::PatchFunction1Types::MappedFile<Type>::MappedFile
fieldTableName_(rhs.fieldTableName_),
pointsName_(rhs.pointsName_),
mapMethod_(rhs.mapMethod_),
filterRadius_(rhs.filterRadius_),
filterSweeps_(rhs.filterSweeps_),
filterFieldPtr_(nullptr),
readerFormat_(rhs.readerFormat_),
readerFile_(rhs.readerFile_),
readerPtr_(nullptr),
mapperPtr_(rhs.mapperPtr_.clone()),
sampleTimes_(rhs.sampleTimes_),
sampleIndex_(rhs.sampleIndex_),
sampleAverage_(rhs.sampleAverage_),
sampleValues_(rhs.sampleValues_),
offset_(rhs.offset_.clone())
{}
{
if (!readerFormat_.empty() && !readerFile_.empty())
{
fileName fName(readerFile_);
fName.expand();
readerPtr_ = surfaceReader::New(readerFormat_, fName);
}
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
@ -165,6 +259,7 @@ void Foam::PatchFunction1Types::MappedFile<Type>::autoMap
}
// Clear interpolator
filterFieldPtr_.reset(nullptr);
mapperPtr_.reset(nullptr);
sampleIndex_ = labelPair(-1, -1);
}
@ -195,6 +290,7 @@ void Foam::PatchFunction1Types::MappedFile<Type>::rmap
}
// Clear interpolator
filterFieldPtr_.reset(nullptr);
mapperPtr_.reset(nullptr);
sampleIndex_ = labelPair(-1, -1);
}
@ -211,13 +307,41 @@ void Foam::PatchFunction1Types::MappedFile<Type>::updateSampledValues
tmp<Field<Type>> tvalues;
// Update sampled data fields
if (readerPtr_)
{
wordList fieldNames = readerPtr_->fieldNames(sampleIndex);
label fieldIndex = fieldNames.find(fieldTableName_);
DebugInfo
<< "checkTable : Update index=" << sampleIndex
<< " field=" << fieldNames[fieldIndex] << endl;
tvalues = readerPtr_->field
(
sampleIndex,
fieldIndex,
pTraits<Type>::zero
);
if (tvalues().size() != mapperPtr_().sourceSize())
{
FatalErrorInFunction
<< "Number of values (" << tvalues().size()
<< ") differs from the number of points ("
<< mapperPtr_().sourceSize() << ")"
<< exit(FatalError);
}
}
else
{
const polyMesh& mesh = this->patch_.boundaryMesh().mesh();
const Time& time = mesh.time();
if (debug)
{
Pout<< "checkTable : Reading values from "
Pout<< "checkTable : Update index=" << sampleIndex
<< " Reading values from "
<<
(
"boundaryData"
@ -269,6 +393,13 @@ void Foam::PatchFunction1Types::MappedFile<Type>::updateSampledValues
tvalues = tmp<Field<Type>>::New(std::move(vals.field()));
}
if (filterFieldPtr_)
{
DebugInfo
<< "apply " << filterSweeps_ << " filter sweeps" << endl;
tvalues = filterFieldPtr_().evaluate(tvalues, filterSweeps_);
}
// From input values to interpolated (sampled) positions
field = mapperPtr_().interpolate(tvalues);
@ -285,7 +416,74 @@ void Foam::PatchFunction1Types::MappedFile<Type>::checkTable
const Time& time = mesh.time();
// Initialise
if (!mapperPtr_)
if (!mapperPtr_ && readerPtr_)
{
auto& reader = readerPtr_();
const meshedSurface& geom = reader.geometry(0);
sampleTimes_ = reader.times();
// We may have been passed in geometry without any faces
// eg, boundaryData
const pointField& samplePoints =
(
geom.nFaces() ? geom.faceCentres() : geom.points()
);
DebugInfo
<< "Read " << samplePoints.size() << " sample points from "
<< readerFile_ << endl
<< "Found times "
<< pointToPointPlanarInterpolation::timeNames(sampleTimes_) << nl;
// tbd: run-time selection
const bool nearestOnly =
(
!mapMethod_.empty() && !mapMethod_.starts_with("planar")
);
// Allocate the interpolator
if (this->faceValues())
{
mapperPtr_.reset
(
new pointToPointPlanarInterpolation
(
samplePoints,
this->localPosition(this->patch_.faceCentres()),
perturb_,
nearestOnly
)
);
}
else
{
mapperPtr_.reset
(
new pointToPointPlanarInterpolation
(
samplePoints,
this->localPosition(this->patch_.localPoints()),
perturb_,
nearestOnly
)
);
}
// Setup median filter (if any)
if (filterSweeps_ > 0)
{
filterFieldPtr_.reset(new FilterField(geom, filterRadius_));
}
else
{
filterFieldPtr_.reset(nullptr);
}
}
else if (!mapperPtr_)
{
// Reread values and interpolate
const fileName samplePointsFile
@ -349,7 +547,6 @@ void Foam::PatchFunction1Types::MappedFile<Type>::checkTable
);
}
// Read the times for which data is available
const fileName samplePointsDir = samplePointsFile.path();
sampleTimes_ = Time::findTimes(samplePointsDir);
@ -358,11 +555,21 @@ void Foam::PatchFunction1Types::MappedFile<Type>::checkTable
<< "Found times "
<< pointToPointPlanarInterpolation::timeNames(sampleTimes_)
<< endl;
// Setup median filter (if any)
if (filterSweeps_ > 0)
{
filterFieldPtr_.reset(new FilterField(samplePoints, filterRadius_));
}
else
{
filterFieldPtr_.reset(nullptr);
}
}
// Find range of current time indices in sampleTimes
Pair<label> timeIndices = instant::findRange
labelPair timeIndices = instant::findRange
(
sampleTimes_,
t, //mesh.time().value(),
@ -555,6 +762,12 @@ void Foam::PatchFunction1Types::MappedFile<Type>::writeEntries
Ostream& os
) const
{
if (!readerFormat_.empty() && !readerFile_.empty())
{
os.writeEntry("readerFormat", readerFormat_);
os.writeEntry("readerFile", readerFile_);
}
os.writeEntryIfDifferent
(
"fieldTable",
@ -579,6 +792,12 @@ void Foam::PatchFunction1Types::MappedFile<Type>::writeEntries
os.writeEntryIfDifferent<scalar>("perturb", 1e-5, perturb_);
if (filterSweeps_ >= 1)
{
os.writeEntry("filterRadius", filterRadius_);
os.writeEntry("filterSweeps", filterSweeps_);
}
if (offset_)
{
offset_->writeData(os);

42
src/meshTools/PatchFunction1/MappedFile/MappedFile.H
View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2018-2021 OpenCFD Ltd.
Copyright (C) 2018-2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -49,10 +49,26 @@ Description
\endtable
Options for the \c mapMethod entry:
\verbatim
\table
Property | Description
nearest | Use nearest points only (avoids triangulation)
planar | Interpolation using 2D Delaunay triangulation
\endverbatim
\endtable
Options for reading and filtering:
\table
Property | Description
sampleFormat | The surfaceReader format (eg, ensight)
sampleFile | <case>/foo/bar/window.case
filterRadius | Search radius [m] for median filter neighbours
filterSweeps | Filter sweeps for median filter
\endtable
Note
The MappedFile handling currently has two different reading types:
a builtin boundaryData reader (default) and a surfaceReader
(specified by sampleFormat and sampleFile keywords). The generic
surfaceReader handling does \em not support field averaging.
SourceFiles
MappedFile.C
@ -66,6 +82,8 @@ SourceFiles
#include "Function1.H"
#include "Pair.H"
#include "instantList.H"
#include "surfaceReader.H"
#include "MappedFileFilterField.H"
#include "pointToPointPlanarInterpolation.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -104,6 +122,24 @@ class MappedFile
//- Interpolation scheme to use (default is empty == "planar")
word mapMethod_;
//- Radius for filter
scalar filterRadius_;
//- Number of median filter sweeps
label filterSweeps_;
//- Median filtering (for input values)
mutable autoPtr<FilterField> filterFieldPtr_;
//- Format name for surfaceReader
word readerFormat_;
//- File name associated with surfaceReader
fileName readerFile_;
//- Meshed surface with fields (for input values)
mutable autoPtr<surfaceReader> readerPtr_;
//- 2D interpolation (for 'planar' mapMethod)
mutable autoPtr<pointToPointPlanarInterpolation> mapperPtr_;

381
src/meshTools/PatchFunction1/MappedFile/MappedFileFilterField.C Normal file
View File

@ -0,0 +1,381 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM 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.
OpenFOAM 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 OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "MappedFileFilterField.H"
#include "oneField.H"
#include "MeshedSurfaces.H"
#include "MinMax.H"
#include "indexedOctree.H"
#include "treeDataPoint.H"
#include <algorithm>
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
namespace PatchFunction1Types
{
defineTypeNameAndDebug(FilterField, 0);
}
}
const Foam::Enum
<
Foam::PatchFunction1Types::FilterField::RBF_type
>
Foam::PatchFunction1Types::FilterField::RBF_typeNames_
{
{ RBF_type::RBF_linear, "linear" },
{ RBF_type::RBF_quadratic, "quadratic" },
{ RBF_type::RBF_linear, "default" },
};
// * * * * * * * * * * * * * * * Local Functions * * * * * * * * * * * * * * //
namespace Foam
{
//- Linear basis function: 1 - (r/r0)
static inline scalar linearWeight
(
const point& p,
const point& p0,
const scalar radiusSqr
)
{
return (1 - ::sqrt(p.distSqr(p0) / radiusSqr));
}
//- Quadratic basis function: 1 - (r/r0)^2
static inline scalar quadraticWeight
(
const point& p,
const point& p0,
const scalar radiusSqr
)
{
return (1 - p.distSqr(p0) / radiusSqr);
}
//- Construct search tree for points
static autoPtr<indexedOctree<treeDataPoint>>
createTree(const pointField& points)
{
// Avoid bounding box error in case of 2D cases
treeBoundBox bb(points);
bb.grow(1e-4);
const int debug = PatchFunction1Types::FilterField::debug;
const int oldDebug = indexedOctreeBase::debug;
if (debug & 2)
{
indexedOctreeBase::debug = 1;
}
auto treePtr = autoPtr<indexedOctree<treeDataPoint>>::New
(
treeDataPoint(points),
bb, // overall search domain
points.size(), // maxLevel
16, // leafsize
1.0 // duplicity
);
indexedOctreeBase::debug = oldDebug;
if (debug & 2)
{
const auto& tree = treePtr();
OFstream os("indexedOctree.obj");
tree.writeOBJ(os);
Info<< "=================" << endl;
Info<< "have " << tree.nodes().size() << " nodes" << nl;
Info<< "=================" << endl;
}
return treePtr;
}
} // End namespace Foam
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template
<
class TreeType,
class RadiusSqrOp,
class BasisFunctionOp,
class PointWeightFieldType
>
void Foam::PatchFunction1Types::FilterField::buildWeightsImpl
(
const TreeType& tree,
const RadiusSqrOp& radiusSqrOp,
const BasisFunctionOp& basisFuncOp,
const PointWeightFieldType& pointWeightFld
)
{
tmp<pointField> tpoints = tree.shapes().centres();
const auto& points = tpoints();
const label nPoints = points.size();
weights_.clear();
addressing_.clear();
weights_.resize(nPoints);
addressing_.resize(nPoints);
labelHashSet neighbours(2*128);
// Catch degenerate case where weight/addressing not actually needed
bool usesNeighbours = false;
for (label pointi = 0; pointi < nPoints; ++pointi)
{
const point& p0 = points[pointi];
auto& currAddr = addressing_[pointi];
auto& currWeight = weights_[pointi];
// Search with local sphere
const scalar radiusSqr = radiusSqrOp(pointi);
tree.findSphere(p0, radiusSqr, neighbours);
// Paranoid, enforce identity weighting as a minimum
if (neighbours.size() < 2)
{
currAddr.resize(1, pointi);
currWeight.resize(1, scalar(1));
continue;
}
usesNeighbours = true;
currAddr = neighbours.sortedToc();
currWeight.resize_nocopy(currAddr.size());
scalar sumWeight = 0;
forAll(currAddr, i)
{
const point& p = points[currAddr[i]];
currWeight[i] = basisFuncOp(p, p0, radiusSqr);
// Eg, face area weighting.
// - cast required for oneField()
currWeight[i] *= static_cast<scalar>(pointWeightFld[i]);
sumWeight += currWeight[i];
}
for (auto& w : currWeight)
{
w /= sumWeight;
}
}
if (!usesNeighbours)
{
// No addressing/weighting required
reset();
}
if (debug && !addressing_.empty())
{
labelMinMax limits(addressing_[0].size());
label total = 0;
for (const labelList& addr : addressing_)
{
const label n = addr.size();
limits.add(n);
total += n;
}
Info<< "Weight neighbours: min=" << limits.min()
<< " avg=" << (total / addressing_.size())
<< " max=" << limits.max() << endl;
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::PatchFunction1Types::FilterField::FilterField
(
const pointField& points,
const scalar radius,
const RBF_type interp
)
{
reset(points, radius, interp);
}
Foam::PatchFunction1Types::FilterField::FilterField
(
const meshedSurface& geom,
const scalar radius,
const bool relative,
const RBF_type interp
)
{
reset(geom, radius, relative, interp);
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::PatchFunction1Types::FilterField::reset()
{
addressing_.clear();
weights_.clear();
}
void Foam::PatchFunction1Types::FilterField::reset
(
const pointField& points,
const scalar radius,
const RBF_type interp
)
{
scalarField dummyWeights;
if (debug)
{
Info<< "Apply " << RBF_typeNames_[interp] << " filter,"
<< " radius=" << radius << nl
<< "Create tree..." << endl;
}
autoPtr<indexedOctree<treeDataPoint>> treePtr
= createTree(points);
const scalar radiusSqr = sqr(radius);
auto weightFunc = linearWeight;
if (interp == RBF_type::RBF_quadratic)
{
weightFunc = quadraticWeight;
}
// Use (constant) absolute size
buildWeightsImpl
(
treePtr(),
[=](label index) -> scalar { return radiusSqr; },
weightFunc,
oneField()
);
}
void Foam::PatchFunction1Types::FilterField::reset
(
const meshedSurface& geom,
const scalar radius,
const bool relative,
const RBF_type interp
)
{
if (!geom.nFaces())
{
// Received geometry without any faces eg, boundaryData
reset(geom.points(), radius, interp);
return;
}
if (debug)
{
Info<< "Apply " << RBF_typeNames_[interp] << " filter,";
if (relative)
{
Info<< " relative";
}
Info<< " radius=" << radius << endl;
Info<< "Create tree..." << endl;
}
autoPtr<indexedOctree<treeDataPoint>> treePtr
= createTree(geom.faceCentres());
const scalar radiusSqr = sqr(radius);
auto weightFunc = linearWeight;
if (interp == RBF_type::RBF_quadratic)
{
weightFunc = quadraticWeight;
}
if (relative)
{
// Use relative size
buildWeightsImpl
(
treePtr(),
[&](label index) -> scalar
{
return (radiusSqr * geom.sphere(index));
},
weightFunc,
geom.magSf()
);
}
else
{
// Use (constant) absolute size
buildWeightsImpl
(
treePtr(),
[=](label index) -> scalar { return radiusSqr; },
weightFunc,
geom.magSf()
);
}
}
// ************************************************************************* //

213
src/meshTools/PatchFunction1/MappedFile/MappedFileFilterField.H Normal file
View File

@ -0,0 +1,213 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM 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.
OpenFOAM 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 OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::PatchFunction1Types::FilterField
Description
The FilterField helper class provides a multi-sweep median filter
for a Field of data associated with a geometric point cloud.
The points can be freestanding or the faceCentres (or points)
of a meshedSurface, for example.
Using an initial specified search radius, the nearest point
neighbours are gathered and addressing/weights are built for them.
This currently uses an area-weighted, linear RBF interpolator
with provision for quadratic RBF interpolator etc.
After the weights and addressing are established,
the evaluate() method can be called to apply a median filter
to data fields, with a specified number of sweeps.
Note
When handling large search radii and/or an extensive number of
filter sweeps, compiling with openmp can yield some speedup.
SourceFiles
MappedFileFilterField.C
MappedFileFilterFieldTemplates.C
\*---------------------------------------------------------------------------*/
#ifndef Foam_PatchFunction1Types_MappedFileFilterField_H
#define Foam_PatchFunction1Types_MappedFileFilterField_H
#include "primitiveFields.H"
#include "pointField.H"
#include "Enum.H"
#include "className.H"
#include "MeshedSurfacesFwd.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace PatchFunction1Types
{
/*---------------------------------------------------------------------------*\
Class FilterField Declaration
\*---------------------------------------------------------------------------*/
class FilterField
{
public:
// Data Types
//- Basis type
enum RBF_type
{
RBF_linear,
RBF_quadratic
};
private:
// Private Data
//- Names for RBF basis types
static const Enum<RBF_type> RBF_typeNames_;
//- Addressing
List<labelList> addressing_;
//- Weights
List<scalarField> weights_;
// Private Member Functions
//- Construct neighbour point weights and addressing based
//- on specified search radius
//
// Search with local sphere
//
// RadiusSqrOp (index) -> radiusSqr
// BasisFunctionOp(point, point, radius^2) -> RBF
// PointWeightFieldType [index] -> weighting
template
<
class TreeType,
class RadiusSqrOp,
class BasisFunctionOp,
class PointWeightFieldType
>
void buildWeightsImpl
(
const TreeType& tree,
const RadiusSqrOp& radiusSqrOp,
const BasisFunctionOp& basisFuncOp,
const PointWeightFieldType& pointWeightFld
);
public:
//- Runtime type information
ClassName("filterField");
// Constructors
//- Default construct
FilterField() noexcept = default;
//- Construct with weights for a field of points
//- (constant search radius)
FilterField
(
const pointField& points,
const scalar radius,
const RBF_type interp = RBF_type::RBF_linear
);
//- Construct with weights for faceCentres of a meshedSurface
//- using a constant search radius search
//- (or relative radius multiplier for the face bounding spheres).
FilterField
(
const meshedSurface& geom,
const scalar radius,
const bool relative = false,
const RBF_type interp = RBF_type::RBF_linear
);
// Member Functions
//- Reset to unweighted (pass-through)
void reset();
//- Create weights for field of points (constant search radius)
void reset
(
const pointField& points,
const scalar radius,
const RBF_type interp = RBF_type::RBF_linear
);
//- Create weights for meshedSurface using a constant search radius
//- or optionally with a search radius multiplier for the
//- face bounding spheres.
void reset
(
const meshedSurface& geom,
const scalar radius,
const bool relative = false,
const RBF_type interp = RBF_type::RBF_linear
);
// Evaluation
//- Return the median smoothed field
template<class Type>
tmp<Field<Type>> evaluate
(
const tmp<Field<Type>>& tinput,
const label nSweeps
) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace PatchFunction1Types
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "MappedFileFilterFieldTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

132
src/meshTools/PatchFunction1/MappedFile/MappedFileFilterFieldTemplates.C Normal file
View File

@ -0,0 +1,132 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM 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.
OpenFOAM 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 OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "indexedOctree.H"
#include "treeDataPoint.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
Foam::tmp<Foam::Field<Type>>
Foam::PatchFunction1Types::FilterField::evaluate
(
const tmp<Field<Type>>& tinput,
const label nSweeps
) const
{
if (nSweeps < 1 || !tinput.good())
{
// Nothing to do
return tinput;
}
label nPoints = tinput().size();
const label nAddr = addressing_.size();
if (!nPoints || !nAddr)
{
// No input or using an identity mapping
return tinput;
}
// Output
auto toutput = tmp<Field<Type>>::New(nPoints);
if (nAddr < nPoints)
{
WarningInFunction
<< "Addressing/weights shorter than input field"
<< endl;
// Restrict addressing space within loop
nPoints = nAddr;
// Fill trailing portion with the input values
toutput.ref().slice(nAddr) = tinput().slice(nAddr);
}
// Intermediate buffer
tmp<Field<Type>> tbuffer;
if (nSweeps == 1)
{
// Simple reference ie enough
tbuffer.cref(tinput);
}
else
{
// Need buffer for swap - copy or steal contents
tbuffer.reset(tinput.ptr());
}
tinput.clear();
// If there are any 'senseless' sweeps
// (ie, with identity mapping of values), they will have
// been detected during addressing construction, can ignore here
for (label sweep = 0; sweep < nSweeps; ++sweep)
{
if (sweep > 0)
{
// Reuse previous output for subsequent sweeps
tbuffer.swap(toutput);
}
const auto& input = tbuffer();
auto& output = toutput.ref();
#pragma omp parallel for if (nPoints > 1000)
for (label pointi = 0; pointi < nPoints; ++pointi)
{
const auto& addr = addressing_[pointi];
const auto& weight = weights_[pointi];
auto& interp = output[pointi];
if (addr.empty())
{
// Could happen if the search radius was really small
interp = input[pointi];
}
else
{
interp = Zero;
forAll(addr, i)
{
interp += (weight[i] * input[addr[i]]);
}
}
}
}
return toutput;
}
// ************************************************************************* //

1
src/surfMesh/Make/files
View File

@ -60,6 +60,7 @@ triSurface/patches/surfacePatch.C
readers = readers
$(readers)/common/surfaceReader.C
$(readers)/boundary/boundaryDataSurfaceReader.C
$(readers)/ensight/ensightSurfaceReader.C
writers = writers

281
src/surfMesh/readers/boundary/boundaryDataSurfaceReader.C Normal file
View File

@ -0,0 +1,281 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM 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.
OpenFOAM 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 OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "boundaryDataSurfaceReader.H"
#include "rawIOField.H"
#include "argList.H"
#include "Time.H"
#include "addToRunTimeSelectionTable.H"
#include <algorithm>
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(boundaryDataSurfaceReader, 0);
addToRunTimeSelectionTable
(
surfaceReader,
boundaryDataSurfaceReader,
fileName
);
}
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
Foam::pointField Foam::boundaryDataSurfaceReader::readPoints
(
const Time& runTime,
const fileName& baseDir,
const word& pointsName
)
{
fileName pointsFile
(
baseDir / (pointsName.empty() ? "points" : pointsName)
);
pointsFile.toAbsolute();
IOobject io
(
pointsFile, // absolute path
runTime,
IOobject::MUST_READ,
IOobject::NO_WRITE,
IOobject::NO_REGISTER,
true // global object (currently not used)
);
DebugInfo
<< "File: " << io.objectPath() << endl;
// Read data (no average value!)
rawIOField<point> rawData(io);
pointField points(std::move(rawData.field()));
DebugInfo
<< "File: " << io.objectPath()
<< " " << points.size() << " points" << endl;
return points;
}
Foam::pointField Foam::boundaryDataSurfaceReader::readPoints
(
const fileName& dirName,
const word& pointsName
)
{
refPtr<Time> timePtr(Time::New(argList::envGlobalPath()));
return readPoints(*timePtr, dirName, pointsName);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::boundaryDataSurfaceReader::readCase()
{
DebugInFunction << endl;
timeValues_ = TimePaths::findTimes(baseDir_);
}
void Foam::boundaryDataSurfaceReader::readGeometry
(
meshedSurface& surf,
const label timeIndex
)
{
surf.clear();
pointField points(this->readPoints(baseDir_, pointsName_));
surf = meshedSurface(std::move(points), faceList());
}
template<class Type>
Foam::tmp<Foam::Field<Type>>
Foam::boundaryDataSurfaceReader::readFieldTemplate
(
const label timeIndex,
const label fieldIndex
) const
{
Type dummyAvg;
return readField<Type>
(
baseDir_,
timeValues_[timeIndex],
fieldNames_[fieldIndex],
dummyAvg
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::boundaryDataSurfaceReader::boundaryDataSurfaceReader
(
const fileName& fName,
const word& pointsName
)
:
surfaceReader(fName),
baseDir_(fName.path()),
pointsName_(pointsName),
timeValues_(),
fieldNames_(),
surfPtr_(nullptr)
{
baseDir_.toAbsolute();
debug = 1;
DebugInFunction << endl;
Info<< "create with " << baseDir_ << endl;
readCase();
}
// * * * * * * * * * * * * * Public Member Functions * * * * * * * * * * * //
const Foam::meshedSurface& Foam::boundaryDataSurfaceReader::geometry
(
const label timeIndex
)
{
DebugInFunction << endl;
if (!surfPtr_)
{
surfPtr_.reset(new meshedSurface);
readGeometry(*surfPtr_, timeIndex);
}
return *surfPtr_;
}
Foam::instantList Foam::boundaryDataSurfaceReader::times() const
{
return timeValues_;
}
Foam::wordList Foam::boundaryDataSurfaceReader::fieldNames
(
const label timeIndex
) const
{
if (timeValues_.empty() || timeIndex >= timeValues_.size())
{
fieldNames_.clear();
return wordList();
}
fileNameList items =
fileHandler().readDir(baseDir_/timeValues_[timeIndex].name());
fieldNames_.resize_nocopy(items.size());
std::transform
(
items.begin(),
items.end(),
fieldNames_.begin(),
[](const fileName& f) { return word(f); }
);
Foam::sort(fieldNames_);
return fieldNames_;
}
Foam::tmp<Foam::Field<Foam::scalar>> Foam::boundaryDataSurfaceReader::field
(
const label timeIndex,
const label fieldIndex,
const scalar& refValue
) const
{
return readFieldTemplate<scalar>(timeIndex, fieldIndex);
}
Foam::tmp<Foam::Field<Foam::vector>> Foam::boundaryDataSurfaceReader::field
(
const label timeIndex,
const label fieldIndex,
const vector& refValue
) const
{
return readFieldTemplate<vector>(timeIndex, fieldIndex);
}
Foam::tmp<Foam::Field<Foam::sphericalTensor>>
Foam::boundaryDataSurfaceReader::field
(
const label timeIndex,
const label fieldIndex,
const sphericalTensor& refValue
) const
{
return readFieldTemplate<sphericalTensor>(timeIndex, fieldIndex);
}
Foam::tmp<Foam::Field<Foam::symmTensor>> Foam::boundaryDataSurfaceReader::field
(
const label timeIndex,
const label fieldIndex,
const symmTensor& refValue
) const
{
return readFieldTemplate<symmTensor>(timeIndex, fieldIndex);
}
Foam::tmp<Foam::Field<Foam::tensor>> Foam::boundaryDataSurfaceReader::field
(
const label timeIndex,
const label fieldIndex,
const tensor& refValue
) const
{
return readFieldTemplate<tensor>(timeIndex, fieldIndex);
}
// ************************************************************************* //

227
src/surfMesh/readers/boundary/boundaryDataSurfaceReader.H Normal file
View File

@ -0,0 +1,227 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM 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.
OpenFOAM 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 OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::boundaryDataSurfaceReader
Description
A boundaryData format surface reader.
However, the "surface" represented by boundaryData is actually only
point data!
// Points
<case>/constant/region0/"boundaryData"/patchName/points
// Values
<case>/constant/region0/"boundaryData"/patchName/TIME/field
SourceFiles
boundaryDataSurfaceReader.C
\*---------------------------------------------------------------------------*/
#ifndef Foam_boundaryDataSurfaceReader_H
#define Foam_boundaryDataSurfaceReader_H
#include "surfaceReader.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// Forward Declarations
class Time;
/*---------------------------------------------------------------------------*\
Class boundaryDataSurfaceReader Declaration
\*---------------------------------------------------------------------------*/
class boundaryDataSurfaceReader
:
public surfaceReader
{
// Private Data
//- Base directory
fileName baseDir_;
//- Name of the 'points' file
word pointsName_;
//- Times
instantList timeValues_;
//- Field names
mutable List<word> fieldNames_;
//- Pointer to the surface
autoPtr<meshedSurface> surfPtr_;
// Private Member Functions
//- Populate time instances
void readCase();
//- Read geometry
void readGeometry(meshedSurface& surf, const label timeIndex);
//- Read and return given field
template<class Type>
tmp<Field<Type>> readFieldTemplate
(
const label timeIndex,
const label fieldIndex
) const;
public:
//- Runtime type information
TypeName("boundaryData");
// Constructors
//- Construct from fileName
explicit boundaryDataSurfaceReader
(
const fileName& fName,
const word& pointsName = "points"
);
//- Destructor
virtual ~boundaryDataSurfaceReader() = default;
// Static Member Functions
//- Read points file
static pointField readPoints
(
const Time& runTime,
const fileName& baseDir,
const word& pointsName = "points"
);
//- Read points file
static pointField readPoints
(
const fileName& baseDir,
const word& pointsName = "points"
);
//- Read and return given field
template<class Type>
static tmp<Field<Type>> readField
(
const Time& runTime,
const fileName& baseDir,
const instant& timeDir,
const word& fieldName,
Type& avg
);
//- Read and return given field
template<class Type>
static tmp<Field<Type>> readField
(
const fileName& baseDir,
const instant& timeDir,
const word& fieldName,
Type& avg
);
// Member Functions
//- Return a reference to the surface geometry
virtual const meshedSurface& geometry(const label timeIndex);
//- Return a list of the available times
virtual instantList times() const;
//- Return a list of the available fields at a given time
virtual wordList fieldNames(const label timeIndex) const;
//- Return a scalar field at a given time
virtual tmp<Field<scalar>> field
(
const label timeIndex,
const label fieldIndex,
const scalar& refValue = pTraits<scalar>::zero
) const;
//- Return a vector field at a given time
virtual tmp<Field<vector>> field
(
const label timeIndex,
const label fieldIndex,
const vector& refValue = pTraits<vector>::zero
) const;
//- Return a sphericalTensor field at a given time
virtual tmp<Field<sphericalTensor>> field
(
const label timeIndex,
const label fieldIndex,
const sphericalTensor& refValue = pTraits<sphericalTensor>::zero
) const;
//- Return a symmTensor field at a given time
virtual tmp<Field<symmTensor>> field
(
const label timeIndex,
const label fieldIndex,
const symmTensor& refValue = pTraits<symmTensor>::zero
) const;
//- Return a tensor field at a given time
virtual tmp<Field<tensor>> field
(
const label timeIndex,
const label fieldIndex,
const tensor& refValue = pTraits<tensor>::zero
) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "boundaryDataSurfaceReaderTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

96
src/surfMesh/readers/boundary/boundaryDataSurfaceReaderTemplates.C Normal file
View File

@ -0,0 +1,96 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM 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.
OpenFOAM 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 OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "boundaryDataSurfaceReader.H"
#include "rawIOField.H"
#include "argList.H"
#include "Time.H"
#include <algorithm>
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
template<class Type>
Foam::tmp<Foam::Field<Type>>
Foam::boundaryDataSurfaceReader::readField
(
const Time& runTime,
const fileName& baseDir,
const instant& timeDir,
const word& fieldName,
Type& avg
)
{
// Reread values and interpolate
fileName valuesFile(baseDir/timeDir.name()/fieldName);
valuesFile.toAbsolute();
IOobject io
(
valuesFile, // absolute path
runTime,
IOobject::MUST_READ,
IOobject::NO_WRITE,
IOobject::NO_REGISTER,
true // global object (currently not used)
);
DebugInfo
<< "File: " << io.objectPath() << endl;
// Read data (TDB: setAverage)
rawIOField<Type> rawData(io, IOobjectOption::READ_IF_PRESENT);
if (rawData.hasAverage())
{
avg = rawData.average();
}
DebugInfo
<< "File: " << io.objectPath()
<< " " << rawData.size() << " values" << endl;
return tmp<Field<Type>>::New(std::move(rawData.field()));
}
template<class Type>
Foam::tmp<Foam::Field<Type>>
Foam::boundaryDataSurfaceReader::readField
(
const fileName& baseDir,
const instant& timeDir,
const word& fieldName,
Type& avg
)
{
refPtr<Time> timePtr(Time::New(argList::envGlobalPath()));
return readField<Type>(baseDir, timeDir, fieldName, avg);
}
// ************************************************************************* //