Base class for noise models. More...
Classes | |
struct | octaveBandInfo |
Octave band information. More... | |
struct | planInfo |
FFTW planner information. More... | |
Public Types | |
enum | weightingType { none, dBA, dBB, dBC, dBD } |
Public Member Functions | |
TypeName ("noiseModel") | |
Runtime type information. More... | |
declareRunTimeSelectionTable (autoPtr, noiseModel, dictionary,(const dictionary &dict, const objectRegistry &obr),(dict, obr)) | |
Run time selection table. More... | |
noiseModel (const noiseModel &)=delete | |
No copy construct. More... | |
void | operator= (const noiseModel &)=delete |
No copy assignment. More... | |
noiseModel (const dictionary &dict, const objectRegistry &obr, const word &name, const bool readFields=true) | |
Constructor. More... | |
virtual | ~noiseModel ()=default |
Destructor. More... | |
virtual bool | read (const dictionary &dict) |
Read from dictionary. More... | |
virtual void | calculate ()=0 |
Abstract call to calculate. More... | |
tmp< Foam::scalarField > | PSD (const scalarField &PSDf) const |
PSD [dB/Hz]. More... | |
tmp< scalarField > | SPL (const scalarField &Prms2, const scalar f) const |
SPL [dB]. More... | |
tmp< scalarField > | SPL (const scalarField &Prms2, const scalarField &f) const |
SPL [dB]. More... | |
void | cleanFFTW () |
Clean up the FFTW. More... | |
void | writeWeightings () const |
Helper function to check weightings. More... | |
Public Member Functions inherited from writeFile | |
writeFile (const objectRegistry &obr, const fileName &prefix, const word &name="undefined", const bool writeToFile=true) | |
Construct from objectRegistry, prefix, fileName. More... | |
writeFile (const objectRegistry &obr, const fileName &prefix, const word &name, const dictionary &dict, const bool writeToFile=true) | |
Construct from objectRegistry, prefix, fileName and read options from dictionary. More... | |
writeFile (const writeFile &wf) | |
Construct copy. More... | |
virtual | ~writeFile ()=default |
Destructor. More... | |
virtual OFstream & | file () |
Return access to the file (if only 1) More... | |
virtual bool | writeToFile () const |
Flag to allow writing to file. More... | |
virtual bool | canWriteToFile () const |
Flag to allow writing to the file. More... | |
virtual bool | canResetFile () const |
Flag to allow resetting the file. More... | |
virtual bool | canWriteHeader () const |
Flag to allow writing the header. More... | |
virtual label | charWidth () const |
Return width of character stream output. More... | |
virtual void | writeCommented (Ostream &os, const string &str) const |
Write a commented string to stream. More... | |
virtual void | writeTabbed (Ostream &os, const string &str) const |
Write a tabbed string to stream. More... | |
virtual void | writeHeader (Ostream &os, const string &str) const |
Write a commented header to stream. More... | |
virtual void | writeCurrentTime (Ostream &os) const |
Write the current time to stream. More... | |
virtual void | writeBreak (Ostream &os) const |
Write a break marker to the stream. More... | |
template<class Type > | |
void | writeHeaderValue (Ostream &os, const string &property, const Type &value) const |
Write a (commented) header property and value pair. More... | |
template<class Type > | |
void | writeValue (Ostream &os, const Type &val) const |
Write a given value to stream with the space delimiter. More... | |
Static Public Member Functions | |
static autoPtr< noiseModel > | New (const dictionary &dict, const objectRegistry &obr) |
Selector. More... | |
Static Public Attributes | |
static const Enum< weightingType > | weightingTypeNames_ |
Static Public Attributes inherited from writeFile | |
static label | addChars = 8 |
Additional characters for writing. More... | |
Protected Member Functions | |
scalar | checkUniformTimeStep (const scalarList ×) const |
Check and return uniform time step. More... | |
bool | validateBounds (const scalarList &p) const |
Return true if all pressure data is within min/max bounds. More... | |
fileName | baseFileDir (const label dataseti) const |
Return the base output directory. More... | |
void | writeFileHeader (Ostream &os, const string &x, const string &y, const UList< Tuple2< string, token >> &headerValues=UList< Tuple2< string, token >>::null()) const |
Write output file header. More... | |
void | writeFreqDataToFile (Ostream &os, const scalarField &f, const scalarField &fx) const |
tmp< scalarField > | uniformFrequencies (const scalar deltaT, const bool check) const |
Create a field of equally spaced frequencies for the current set of data - assumes a constant time step. More... | |
tmp< scalarField > | octaves (const scalarField &data, const scalarField &f, const labelUList &freqBandIDs) const |
Generate octave data. More... | |
tmp< scalarField > | Pf (const scalarField &p) const |
Return the fft of the given pressure data. More... | |
tmp< scalarField > | meanPf (const scalarField &p) const |
Return the multi-window mean fft of the complete pressure data [Pa]. More... | |
tmp< scalarField > | RMSmeanPf (const scalarField &p) const |
Return the multi-window RMS mean fft of the complete pressure data [Pa]. More... | |
tmp< scalarField > | PSDf (const scalarField &p, const scalar deltaT) const |
Return the multi-window Power Spectral Density (PSD) of the complete pressure data [Pa^2/Hz]. More... | |
scalar | RAf (const scalar f) const |
A weighting function. More... | |
scalar | gainA (const scalar f) const |
A weighting as gain in dB. More... | |
scalar | RBf (const scalar f) const |
B weighting function. More... | |
scalar | gainB (const scalar f) const |
B weighting as gain in dB. More... | |
scalar | RCf (const scalar f) const |
C weighting function. More... | |
scalar | gainC (const scalar f) const |
C weighting as gain in dB. More... | |
scalar | RDf (const scalar f) const |
D weighting function. More... | |
scalar | gainD (const scalar f) const |
D weighting as gain in dB. More... | |
Protected Member Functions inherited from writeFile | |
void | initStream (Ostream &os) const |
Initialise the output stream for writing. More... | |
fileName | baseFileDir () const |
Return the base directory for output. More... | |
fileName | baseTimeDir () const |
Return the base directory for the current time value. More... | |
fileName | filePath (const fileName &fName) const |
Return the full path for the supplied file name. More... | |
virtual autoPtr< OFstream > | newFile (const fileName &fName) const |
Return autoPtr to a new file using file name. More... | |
virtual autoPtr< OFstream > | newFileAtTime (const word &name, scalar timeValue) const |
Return autoPtr to a new file for a given time. More... | |
virtual autoPtr< OFstream > | newFileAtStartTime (const word &name) const |
Return autoPtr to a new file using the simulation start time. More... | |
virtual void | resetFile (const word &name) |
Reset internal file pointer to new file with new name. More... | |
Omanip< int > | valueWidth (const label offset=0) const |
Return the value width when writing to stream with optional offset. More... | |
void | operator= (const writeFile &)=delete |
No copy assignment. More... | |
virtual autoPtr< OFstream > | createFile (const word &name, scalar timeValue) const |
Deprecated(2022-09) Return autoPtr to a new file for a given time. More... | |
virtual autoPtr< OFstream > | createFile (const word &name) const |
Deprecated(2022-09) Return autoPtr to a new file using the simulation start time. More... | |
Static Protected Member Functions | |
static label | findStartTimeIndex (const instantList &allTimes, const scalar startTime) |
Find and return start time index. More... | |
static void | setOctaveBands (const scalarField &f, const scalar fLower, const scalar fUpper, const scalar octave, labelList &fBandIDs, scalarField &fCentre) |
Return a list of the frequency indices wrt f field that correspond to the bands limits for a given octave. More... | |
Protected Attributes | |
const dictionary | dict_ |
Copy of dictionary used for construction. More... | |
scalar | rhoRef_ |
Reference density (to convert from kinematic to static pressure) More... | |
label | nSamples_ |
Number of samples in sampling window, default = 2^16. More... | |
scalar | fLower_ |
Lower frequency limit, default = 25Hz. More... | |
scalar | fUpper_ |
Upper frequency limit, default = 10kHz. More... | |
scalar | sampleFreq_ |
Prescribed sample frequency. More... | |
scalar | startTime_ |
Start time, default = 0s. More... | |
autoPtr< windowModel > | windowModelPtr_ |
Window model. More... | |
weightingType | SPLweighting_ |
Weighting. More... | |
scalar | dBRef_ |
Reference for dB calculation, default = 2e-5. More... | |
scalar | minPressure_ |
Min pressure value. More... | |
scalar | maxPressure_ |
Min pressure value. More... | |
fileName | outputPrefix_ |
Output file prefix, default = ''. More... | |
bool | writePrmsf_ |
Write Prmsf; default = yes. More... | |
bool | writeSPL_ |
Write SPL; default = yes. More... | |
bool | writePSD_ |
Write PSD; default = yes. More... | |
bool | writePSDf_ |
Write PSDf; default = yes. More... | |
bool | writeOctaves_ |
Write writeOctaves; default = yes. More... | |
planInfo | planInfo_ |
Plan information for FFTW. More... | |
Protected Attributes inherited from writeFile | |
const objectRegistry & | fileObr_ |
Reference to the region objectRegistry. More... | |
const fileName | prefix_ |
Prefix. More... | |
word | fileName_ |
Name of file. More... | |
autoPtr< OFstream > | filePtr_ |
File pointer. More... | |
label | writePrecision_ |
Write precision. More... | |
bool | writeToFile_ |
Flag to enable/disable writing to file. More... | |
bool | updateHeader_ |
Flag to update the header, e.g. on mesh changes. Default is true. More... | |
bool | writtenHeader_ |
Flag to identify whether the header has been written. More... | |
bool | useUserTime_ |
Flag to use the specified user time, e.g. CA deg instead of seconds. Default = true. More... | |
scalar | startTime_ |
Start time value. More... | |
Base class for noise models.
Data is read from a dictionary, e.g.
rhoRef 1; N 4096; minFreq 25; maxFreq 10000; startTime 0; outputPrefix "test1"; SPLweighting dBA; // Optional write options dictionary writeOptions { writePrmsf no; writeSPL yes; writePSD yes; writePSDf no; writeOctaves yes; }
where
Property | Description | Required | Default value |
---|---|---|---|
rhoRef | Reference density | no | 1 |
N | Number of samples in sampling window | no | 65536 (2^16) |
minFreq | Lower frequency bounds (fl) | no | 25 |
maxFreq | Upper frequency bounds (fu) | no | 10000 |
sampleFreq | Sample frequency | no | 0 |
startTime | Start time | no | 0 |
outputPrefix | Prefix applied to output files | no | '' |
SPLweighting | Weighting: dBA, dBB, dBC, DBD | no | none |
dBRef | Reference for dB calculation | no | 2e-5 |
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 |
If provided, the sampleFreq is used to define the deltaT between samples. Otherwise the deltaT is determined from the time samples themselves.
Definition at line 175 of file noiseModel.H.
|
strong |
Enumerator | |
---|---|
none | |
dBA | |
dBB | |
dBC | |
dBD |
Definition at line 181 of file noiseModel.H.
|
delete |
No copy construct.
noiseModel | ( | const dictionary & | dict, |
const objectRegistry & | obr, | ||
const word & | name, | ||
const bool | readFields = true |
||
) |
Constructor.
Definition at line 605 of file noiseModel.C.
References Foam::ensightOutput::debug, dict, Foam::read(), and Foam::readFields().
|
virtualdefault |
Destructor.
|
protected |
Check and return uniform time step.
Definition at line 177 of file noiseModel.C.
References UList< T >::back(), Foam::constant::electromagnetic::e, Foam::endl(), Foam::exit(), Foam::FatalError, FatalErrorInFunction, UList< T >::front(), Foam::mag(), UList< T >::size(), and WarningInFunction.
Referenced by surfaceNoise::initialise().
|
protected |
Return true if all pressure data is within min/max bounds.
Definition at line 216 of file noiseModel.C.
References Foam::endl(), forAll, noiseModel::maxPressure_, noiseModel::minPressure_, Foam::nl, p, and WarningInFunction.
|
inlinestaticprotected |
Find and return start time index.
Definition at line 350 of file noiseModel.H.
References instant::findStart(), and startTime.
|
protected |
Return the base output directory.
Definition at line 238 of file noiseModel.C.
References Foam::name(), and Foam::type().
|
protected |
Write output file header.
Definition at line 250 of file noiseModel.C.
References Foam::nl, os(), Foam::writeHeader(), x, and y.
Referenced by surfaceNoise::calculate().
|
protected |
Definition at line 278 of file noiseModel.C.
References f(), forAll, Foam::nl, os(), and Foam::tab.
Referenced by surfaceNoise::calculate().
|
protected |
Create a field of equally spaced frequencies for the current set of data - assumes a constant time step.
Definition at line 295 of file noiseModel.C.
References Foam::check(), Foam::endl(), f(), forAll, N(), tmp< T >::New(), WarningInFunction, and Foam::Zero.
Referenced by surfaceNoise::calculate().
|
staticprotected |
Return a list of the frequency indices wrt f field that correspond to the bands limits for a given octave.
Definition at line 48 of file noiseModel.C.
References DynamicList< T, SizeMin >::append(), Foam::endl(), f(), Foam::flatOutput(), forAll, Foam::pow(), UList< T >::size(), List< T >::transfer(), and WarningInFunction.
Referenced by surfaceNoise::calculate().
|
protected |
Generate octave data.
Definition at line 332 of file noiseModel.C.
References Foam::endl(), f(), tmp< T >::New(), UList< T >::size(), WarningInFunction, and Foam::Zero.
Referenced by surfaceNoise::calculate().
|
protected |
Return the fft of the given pressure data.
Definition at line 383 of file noiseModel.C.
References Foam::abort(), Foam::FatalError, FatalErrorInFunction, forAll, Foam::mag(), n, tmp< T >::New(), p, Foam::constant::atomic::re, fft::realTransform1D(), UList< T >::size(), and Foam::sqrt().
|
protected |
Return the multi-window mean fft of the complete pressure data [Pa].
Definition at line 426 of file noiseModel.C.
References N(), tmp< T >::New(), p, and Foam::Zero.
|
protected |
Return the multi-window RMS mean fft of the complete pressure data [Pa].
Definition at line 447 of file noiseModel.C.
References N(), p, Foam::sqr(), Foam::sqrt(), and Foam::Zero.
Referenced by surfaceNoise::calculate().
|
protected |
Return the multi-window Power Spectral Density (PSD) of the complete pressure data [Pa^2/Hz].
Definition at line 466 of file noiseModel.C.
References Foam::average(), Foam::ensightOutput::debug, Foam::endl(), N(), tmp< T >::New(), p, Foam::Pout, Foam::sqr(), and Foam::Zero.
Referenced by surfaceNoise::calculate().
|
protected |
A weighting function.
Definition at line 501 of file noiseModel.C.
References Foam::constant::physicoChemical::c1, Foam::constant::physicoChemical::c2, f(), Foam::sqr(), and Foam::sqrt().
|
protected |
A weighting as gain in dB.
Definition at line 518 of file noiseModel.C.
References f(), and Foam::log10().
|
protected |
B weighting function.
Definition at line 529 of file noiseModel.C.
References Foam::constant::physicoChemical::c1, Foam::constant::physicoChemical::c2, f(), Foam::sqr(), and Foam::sqrt().
|
protected |
B weighting as gain in dB.
Definition at line 545 of file noiseModel.C.
References f(), and Foam::log10().
|
protected |
C weighting function.
Definition at line 556 of file noiseModel.C.
References Foam::constant::physicoChemical::c1, Foam::constant::physicoChemical::c2, f(), and Foam::sqr().
|
protected |
C weighting as gain in dB.
Definition at line 567 of file noiseModel.C.
References f(), and Foam::log10().
|
protected |
D weighting function.
Definition at line 578 of file noiseModel.C.
References f(), Foam::sqr(), and Foam::sqrt().
|
protected |
D weighting as gain in dB.
Definition at line 591 of file noiseModel.C.
References f(), and Foam::log10().
TypeName | ( | "noiseModel" | ) |
Runtime type information.
declareRunTimeSelectionTable | ( | autoPtr | , |
noiseModel | , | ||
dictionary | , | ||
(const dictionary &dict, const objectRegistry &obr) | , | ||
(dict, obr) | |||
) |
Run time selection table.
|
delete |
No copy assignment.
|
static |
Selector.
Definition at line 26 of file noiseModelNew.C.
References dict, Foam::endl(), Foam::exit(), Foam::FatalIOError, FatalIOErrorInLookup, and Foam::Info.
|
virtual |
Read from dictionary.
Reimplemented from writeFile.
Reimplemented in surfaceNoise, and pointNoise.
Definition at line 649 of file noiseModel.C.
References dict, Foam::endl(), Foam::exit(), Foam::FatalIOError, FatalIOErrorInFunction, Foam::Info, windowModel::New(), Foam::nl, writeFile::read(), dictionary::readIfPresent(), dictionary::readIfPresentCompat(), Foam::readWriteOption(), and dictionary::subOrEmptyDict().
Referenced by pointNoise::read(), and surfaceNoise::read().
|
pure virtual |
Abstract call to calculate.
Implemented in surfaceNoise, and pointNoise.
Foam::tmp< Foam::scalarField > PSD | ( | const scalarField & | PSDf | ) | const |
PSD [dB/Hz].
Definition at line 758 of file noiseModel.C.
References Foam::safeLog10(), and Foam::sqr().
Referenced by surfaceNoise::calculate().
Foam::tmp< Foam::scalarField > SPL | ( | const scalarField & | Prms2, |
const scalar | f | ||
) | const |
SPL [dB].
Definition at line 767 of file noiseModel.C.
References Foam::abort(), f(), Foam::FatalError, FatalErrorInFunction, tmp< T >::ref(), Foam::safeLog10(), and Foam::sqr().
Referenced by surfaceNoise::calculate().
Foam::tmp< Foam::scalarField > SPL | ( | const scalarField & | Prms2, |
const scalarField & | f | ||
) | const |
SPL [dB].
Definition at line 814 of file noiseModel.C.
References Foam::abort(), f(), Foam::FatalError, FatalErrorInFunction, forAll, tmp< T >::ref(), Foam::safeLog10(), and Foam::sqr().
void cleanFFTW | ( | ) |
Clean up the FFTW.
Definition at line 872 of file noiseModel.C.
void writeWeightings | ( | ) | const |
Helper function to check weightings.
Definition at line 883 of file noiseModel.C.
|
static |
Definition at line 190 of file noiseModel.H.
|
protected |
Copy of dictionary used for construction.
Definition at line 228 of file noiseModel.H.
Referenced by pointNoise::calculate().
|
protected |
Reference density (to convert from kinematic to static pressure)
Definition at line 233 of file noiseModel.H.
|
protected |
Number of samples in sampling window, default = 2^16.
Definition at line 238 of file noiseModel.H.
|
protected |
Lower frequency limit, default = 25Hz.
Definition at line 243 of file noiseModel.H.
Referenced by surfaceNoise::calculate().
|
protected |
Upper frequency limit, default = 10kHz.
Definition at line 248 of file noiseModel.H.
Referenced by surfaceNoise::calculate().
|
protected |
Prescribed sample frequency.
Definition at line 253 of file noiseModel.H.
Referenced by surfaceNoise::initialise().
|
protected |
Start time, default = 0s.
Definition at line 258 of file noiseModel.H.
Referenced by surfaceNoise::initialise().
|
protected |
Window model.
Definition at line 263 of file noiseModel.H.
Referenced by surfaceNoise::calculate(), and surfaceNoise::initialise().
|
protected |
Weighting.
Definition at line 268 of file noiseModel.H.
|
protected |
Reference for dB calculation, default = 2e-5.
Definition at line 273 of file noiseModel.H.
|
protected |
Min pressure value.
Definition at line 281 of file noiseModel.H.
Referenced by noiseModel::validateBounds().
|
protected |
Min pressure value.
Definition at line 286 of file noiseModel.H.
Referenced by noiseModel::validateBounds().
|
protected |
Output file prefix, default = ''.
Definition at line 294 of file noiseModel.H.
|
protected |
Write Prmsf; default = yes.
Definition at line 299 of file noiseModel.H.
Referenced by surfaceNoise::calculate().
|
protected |
Write SPL; default = yes.
Definition at line 304 of file noiseModel.H.
Referenced by surfaceNoise::calculate().
|
protected |
Write PSD; default = yes.
Definition at line 309 of file noiseModel.H.
Referenced by surfaceNoise::calculate().
|
protected |
Write PSDf; default = yes.
Definition at line 314 of file noiseModel.H.
Referenced by surfaceNoise::calculate().
|
protected |
Write writeOctaves; default = yes.
Definition at line 319 of file noiseModel.H.
Referenced by surfaceNoise::calculate().
|
mutableprotected |
Plan information for FFTW.
Definition at line 327 of file noiseModel.H.