Commit de3f38da authored by Andrew Heather's avatar Andrew Heather
Browse files

Merge branch 'feature-noise-additions' into 'develop'

ENH: Feature noise additions

See merge request !98
parents a886410d a7a3126a
......@@ -28,6 +28,7 @@ License
#include "IOmanip.H"
#include "Pair.H"
#include "OSspecific.H"
#include "SubField.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
......@@ -160,6 +161,45 @@ Foam::scalarField& Foam::graph::y()
}
void Foam::graph::setXRange(const scalar x0, const scalar x1)
{
if (x1 < x0)
{
FatalErrorInFunction
<< "When setting limits, x1 must be greater than x0" << nl
<< " x0: " << x0 << nl
<< " x1: " << x1 << nl
<< abort(FatalError);
}
label i0 = 0;
label i1 = 0;
forAll(x_, i)
{
if (x_[i] < x0)
{
i0 = i + 1;
}
if (x_[i] < x1)
{
i1 = i;
}
}
label nX = i1 - i0 + 1;
scalarField xNew(SubField<scalar>(x_, nX, i0));
x_.transfer(xNew);
forAllIter(HashPtrTable<curve>, *this, iter)
{
curve* c = iter();
scalarField cNew(SubField<scalar>(*c, nX, i0));
c->transfer(cNew);
}
}
Foam::autoPtr<Foam::graph::writer> Foam::graph::writer::New
(
const word& graphFormat
......
......@@ -178,6 +178,9 @@ public:
scalarField& y();
// Limit the data for all axes based on x-limits
void setXRange(const scalar x0, const scalar x1);
// Write
......
......@@ -353,7 +353,7 @@ Foam::functionObjects::forces::devRhoReff() const
else if (foundObject<dictionary>("transportProperties"))
{
const dictionary& transportProperties =
lookupObject<dictionary>("transportProperties");
lookupObject<dictionary>("transportProperties");
dimensionedScalar nu(transportProperties.lookup("nu"));
......
......@@ -142,7 +142,7 @@ public:
// Member Functions
//- Are angles in degrees?
bool inDegrees() const;
bool inDegrees() const;
//- Non-const access to inDegrees
bool& inDegrees();
......
......@@ -35,6 +35,26 @@ namespace Foam
// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
void Foam::noiseModel::readWriteOption
(
const dictionary& dict,
const word& lookup,
bool& option
) const
{
dict.readIfPresent(lookup, option);
if (option)
{
Info<< " " << lookup << ": " << "yes" << endl;
}
else
{
Info<< " " << lookup << ": " << "no" << endl;
}
}
Foam::scalar Foam::noiseModel::checkUniformTimeStep
(
const scalarList& times
......@@ -92,6 +112,22 @@ Foam::label Foam::noiseModel::findStartTimeIndex
}
Foam::fileName Foam::noiseModel::baseFileDir(const label dataseti) const
{
fileName baseDir("$FOAM_CASE");
word datasetName("input" + Foam::name(dataseti));
baseDir =
baseDir.expand()
/"postProcessing"
/"noise"
/outputPrefix_
/type()
/datasetName;
return baseDir;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::noiseModel::noiseModel(const dictionary& dict, const bool readFields)
......@@ -101,9 +137,16 @@ Foam::noiseModel::noiseModel(const dictionary& dict, const bool readFields)
nSamples_(65536),
fLower_(25),
fUpper_(10000),
customBounds_(false),
startTime_(0),
windowModelPtr_(),
graphFormat_("raw")
graphFormat_("raw"),
outputPrefix_(),
writePrmsf_(true),
writeSPL_(true),
writePSD_(true),
writePSDf_(true),
writeOctaves_(true)
{
if (readFields)
{
......@@ -124,10 +167,18 @@ bool Foam::noiseModel::read(const dictionary& dict)
{
dict.readIfPresent("rhoRef", rhoRef_);
dict.readIfPresent("N", nSamples_);
dict.readIfPresent("fl", fLower_);
dict.readIfPresent("fu", fUpper_);
customBounds_ = false;
if (dict.readIfPresent("fl", fLower_))
{
customBounds_ = true;
}
if (dict.readIfPresent("fu", fUpper_))
{
customBounds_ = true;
}
dict.readIfPresent("startTime", startTime_);
dict.readIfPresent("graphFormat", graphFormat_);
dict.readIfPresent("outputPrefix", outputPrefix_);
// Check number of samples - must be a power of 2 for our FFT
bool powerOf2 = ((nSamples_ != 0) && !(nSamples_ & (nSamples_ - 1)));
......@@ -164,6 +215,14 @@ bool Foam::noiseModel::read(const dictionary& dict)
}
Info<< " Write options:" << endl;
dictionary optDict(dict.subOrEmptyDict("writeOptions"));
readWriteOption(optDict, "writePrmsf", writePrmsf_);
readWriteOption(optDict, "writeSPL", writeSPL_);
readWriteOption(optDict, "writePSD", writePSD_);
readWriteOption(optDict, "writeOctaves", writeOctaves_);
windowModelPtr_ = windowModel::New(dict, nSamples_);
return true;
......
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2015-2016 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2015-2017 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -35,17 +35,35 @@ Description
fl 25;
fu 25;
startTime 0;
outputPrefix "test1";
// Optional write options dictionary
writeOptions
{
writePrmsf no;
writeSPL yes;
writePSD yes;
writePSDf no;
writeOctaves yes;
}
\endverbatim
where
\table
Property | Description | Required | Default value
rhoRef | Reference density | no | 1
N | Number of samples in sampling window | no | 65536 (2^16)
fl | Lower frequency bounds | no | 25
fu | Upper frequency bounds | no | 10000
startTime | Start time | no | 0
graphFormat | Graph format | no | raw
Property | Description | Required | Default value
rhoRef | Reference density | no | 1
N | Number of samples in sampling window | no | 65536 (2^16)
fl | Lower frequency bounds | no | 25
fu | Upper frequency bounds | no | 10000
startTime | Start time | no | 0
outputPrefix | Prefix applied to output files| no | ''
graphFormat | Graph format | no | raw
writePrmsf | Write Prmsf data | no | yes
writeSPL | Write SPL data | no | yes
writePSD | Write PSD data | no | yes
writePSDf | Write PSDf data | no | yes
writeOctaves | Write octaves data | no | yes
\endtable
Note
......@@ -108,6 +126,9 @@ protected:
//- Upper frequency limit, default = 10kHz
scalar fUpper_;
//- Flag to indicate that custom frequenct bounds are being used
bool customBounds_;
//- Start time, default = 0s
scalar startTime_;
......@@ -118,8 +139,37 @@ protected:
word graphFormat_;
// Write options
//- Output file prefix, default = ''
fileName outputPrefix_;
//- Write Prmsf; default = yes
bool writePrmsf_;
//- Write SPL; default = yes
bool writeSPL_;
//- Write PSD; default = yes
bool writePSD_;
//- Write PSDf; default = yes
bool writePSDf_;
//- Write writeOctaves; default = yes
bool writeOctaves_;
// Protected Member Functions
//- Helper function to read write options and provide info feedback
void readWriteOption
(
const dictionary& dict,
const word& lookup,
bool& option
) const;
//- Check and return uniform time step
scalar checkUniformTimeStep
(
......@@ -133,6 +183,9 @@ protected:
const scalar startTime
) const;
//- Return the base output directory
fileName baseFileDir(const label dataseti) const;
public:
......
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2015-2016 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2015-2017 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -66,7 +66,11 @@ void pointNoise::filterTimeData
}
void pointNoise::processData(const Function1Types::CSV<scalar>& data)
void pointNoise::processData
(
const label dataseti,
const Function1Types::CSV<scalar>& data
)
{
Info<< "Reading data file " << data.fName() << endl;
......@@ -86,7 +90,7 @@ void pointNoise::processData(const Function1Types::CSV<scalar>& data)
windowModelPtr_->validate(t.size());
const windowModel& win = windowModelPtr_();
const scalar deltaf = 1.0/(deltaT*win.nSamples());
fileName outDir(fileName("postProcessing")/"noise"/typeName/fNameBase);
fileName outDir(baseFileDir(dataseti)/fNameBase);
// Create the fft
noiseFFT nfft(deltaT, p);
......@@ -97,13 +101,27 @@ void pointNoise::processData(const Function1Types::CSV<scalar>& data)
// RMS pressure [Pa]
graph Prmsf(nfft.RMSmeanPf(win));
Info<< " Creating graph for " << Prmsf.title() << endl;
Prmsf.write(outDir, graph::wordify(Prmsf.title()), graphFormat_);
if (customBounds_)
{
Prmsf.setXRange(fLower_, fUpper_);
}
if (writePrmsf_)
{
Info<< " Creating graph for " << Prmsf.title() << endl;
Prmsf.write(outDir, graph::wordify(Prmsf.title()), graphFormat_);
}
// PSD [Pa^2/Hz]
graph PSDf(nfft.PSDf(win));
Info<< " Creating graph for " << PSDf.title() << endl;
PSDf.write(outDir, graph::wordify(PSDf.title()), graphFormat_);
if (customBounds_)
{
PSDf.setXRange(fLower_, fUpper_);
}
if (writePSDf_)
{
Info<< " Creating graph for " << PSDf.title() << endl;
PSDf.write(outDir, graph::wordify(PSDf.title()), graphFormat_);
}
// PSD [dB/Hz]
graph PSDg
......@@ -114,8 +132,12 @@ void pointNoise::processData(const Function1Types::CSV<scalar>& data)
Prmsf.x(),
noiseFFT::PSD(PSDf.y())
);
Info<< " Creating graph for " << PSDg.title() << endl;
PSDg.write(outDir, graph::wordify(PSDg.title()), graphFormat_);
if (writePSD_)
{
Info<< " Creating graph for " << PSDg.title() << endl;
PSDg.write(outDir, graph::wordify(PSDg.title()), graphFormat_);
}
// SPL [dB]
graph SPLg
......@@ -126,79 +148,58 @@ void pointNoise::processData(const Function1Types::CSV<scalar>& data)
Prmsf.x(),
noiseFFT::SPL(PSDf.y()*deltaf)
);
Info<< " Creating graph for " << SPLg.title() << endl;
SPLg.write(outDir, graph::wordify(SPLg.title()), graphFormat_);
labelList octave13BandIDs;
scalarField octave13FreqCentre;
noiseFFT::octaveBandInfo
(
Prmsf.x(),
fLower_,
fUpper_,
3,
octave13BandIDs,
octave13FreqCentre
);
// 1/3 octave data
// ---------------
// PSD [Pa^2/Hz]
graph PSD13f(nfft.octaves(PSDf, octave13BandIDs, false));
// Integrated PSD = P(rms)^2 [Pa^2]
graph Prms13f2(nfft.octaves(PSDf, octave13BandIDs, true));
graph PSD13g
(
"PSD13_dB_Hz(fm)",
"fm [Hz]",
"PSD(fm) [dB_Hz]",
octave13FreqCentre,
noiseFFT::PSD(PSD13f.y())
);
Info<< " Creating graph for " << PSD13g.title() << endl;
PSD13g.write(outDir, graph::wordify(PSD13g.title()), graphFormat_);
graph SPL13g
(
"SPL13_dB(fm)",
"fm [Hz]",
"SPL(fm) [dB]",
octave13FreqCentre,
noiseFFT::SPL(Prms13f2.y())
);
Info<< " Creating graph for " << SPL13g.title() << endl;
SPL13g.write(outDir, graph::wordify(SPL13g.title()), graphFormat_);
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void pointNoise::calculate()
{
// Point data only handled by master
if (!Pstream::master())
if (writeSPL_)
{
return;
Info<< " Creating graph for " << SPLg.title() << endl;
SPLg.write(outDir, graph::wordify(SPLg.title()), graphFormat_);
}
if (inputFileNames_.size())
{
forAll(inputFileNames_, i)
{
const fileName fName = inputFileNames_[i].expand();
Function1Types::CSV<scalar> data("pressure", dict_, "Data", fName);
processData(data);
}
}
else
if (writeOctaves_)
{
Function1Types::CSV<scalar> data("pressure", dict_, "Data");
processData(data);
labelList octave13BandIDs;
scalarField octave13FreqCentre;
noiseFFT::octaveBandInfo
(
Prmsf.x(),
fLower_,
fUpper_,
3,
octave13BandIDs,
octave13FreqCentre
);
// 1/3 octave data
// ---------------
// PSD [Pa^2/Hz]
graph PSD13f(nfft.octaves(PSDf, octave13BandIDs, false));
// Integrated PSD = P(rms)^2 [Pa^2]
graph Prms13f2(nfft.octaves(PSDf, octave13BandIDs, true));
graph PSD13g
(
"PSD13_dB_Hz(fm)",
"fm [Hz]",
"PSD(fm) [dB_Hz]",
octave13FreqCentre,
noiseFFT::PSD(PSD13f.y())
);
Info<< " Creating graph for " << PSD13g.title() << endl;
PSD13g.write(outDir, graph::wordify(PSD13g.title()), graphFormat_);
graph SPL13g
(
"SPL13_dB(fm)",
"fm [Hz]",
"SPL(fm) [dB]",
octave13FreqCentre,
noiseFFT::SPL(Prms13f2.y())
);
Info<< " Creating graph for " << SPL13g.title() << endl;
SPL13g.write(outDir, graph::wordify(SPL13g.title()), graphFormat_);
}
}
......@@ -222,11 +223,45 @@ pointNoise::~pointNoise()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void pointNoise::calculate()
{
// Point data only handled by master
if (!Pstream::master())
{
return;
}
forAll(inputFileNames_, filei)
{
fileName fName = inputFileNames_[filei];
fName.expand();
if (!fName.isAbsolute())
{
fName = "$FOAM_CASE"/fName;
}
fName.expand();
Function1Types::CSV<scalar> data("pressure", dict_, "Data", fName);
processData(filei, data);
}
}
bool pointNoise::read(const dictionary& dict)
{
if (noiseModel::read(dict))
{
dict.readIfPresent("inputFiles", inputFileNames_);
if (!dict.readIfPresent("inputFiles", inputFileNames_))
{
inputFileNames_.setSize(1);
// Note: unsafe lookup of file name from pressureData sub-dict!
dict.subDict("pressureData").lookup("fileName")
>> inputFileNames_[0];
}
return true;
}
......
......@@ -117,7 +117,11 @@ protected:
) const;
//- Process the CSV data
void processData(const Function1Types::CSV<scalar>& data);
void processData
(
const label dataseti,
const Function1Types::CSV<scalar>& data
);
public:
......
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2015-2016 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2015-2017 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -226,20 +226,18 @@ void surfaceNoise::readSurfaceData
Foam::scalar surfaceNoise::writeSurfaceData
(
const fileName& outDirBase,
const word& fName,
const word& groupName,
const word& title,
const scalar freq,
const scalarField& data,
const labelList& procFaceOffset
const labelList& procFaceOffset,
const bool writeSurface
) const
{
Info<< " processing " << title << " for frequency " << freq << endl;
fileName outDir
(
fileName("postProcessing")/"noise"/groupName/Foam::name(freq)
);
const fileName outDir(outDirBase/Foam::name(freq));
if (Pstream::parRun())
{
......@@ -279,20 +277,23 @@ Foam::scalar surfaceNoise::writeSurfaceData
}
}
// could also have meshedSurface implement meshedSurf
fileName outFileName = writerPtr_->write
(
outDir,
fName,
meshedSurfRef
// Could also have meshedSurface implement meshedSurf
if (writeSurface)
{
fileName outFileName = writerPtr_->write
(