Calculate Rayleigh ellipticity curves. More...
#include <Ellipticity.h>
Public Member Functions | |
Ellipticity (int nModes, const QVector< double > *x) | |
double | peakMisfit (int modeIndex, const RealStatisticalValue &val, Dispersion &disp, Rayleigh *waveModel) |
QList< double > | peaks (int modeIndex, Dispersion &disp, Rayleigh *waveModel) |
Calculate Rayleigh ellipticity curves.
The main function used to calculate the ellipticity curve is calculate().
Experience has shown that a 1e-10 precision on dc is necessary to refine the ellipticity down to 1e-3 The precision used to compute the original curve must be the same as for the refines.
QGpCoreWave::Ellipticity::Ellipticity | ( | int | nModes, |
const QVector< double > * | x | ||
) | [inline] |
: ModalStorage(nModes, x) {}
double QGpCoreWave::Ellipticity::peakMisfit | ( | int | modeIndex, |
const RealStatisticalValue & | val, | ||
Dispersion & | disp, | ||
Rayleigh * | waveModel | ||
) |
Find the frequency of the closest peak to val and return the misfit
References QGpCoreTools::endl(), QGpCoreTools::Value< numberType >::isValid(), QGpCoreTools::StatisticalValue< numberType >::mean(), misfit(), QGpCoreWave::ModalStorage::mode(), QGpCoreWave::RayleighTemplate< RealType >::model(), QGpCoreWave::ModalStorage::refineSort(), QGpCoreWave::Seismic1DModel::roughFundamentalFrequency(), QGpCoreTools::StatisticalValue< numberType >::stddev(), QGpCoreTools::tr(), TRACE, QGpCoreWave::ModalStorage::x(), and QGpCoreWave::ModalStorage::xCount().
Referenced by QGpCoreWave::EllipticityFactory::peakMisfit().
{ TRACE; double f0, devf0, misfit=1e99; f0=val.mean()*2.0*M_PI; if(val.stddev()>0) devf0=1.0/ (val.stddev()*2.0*M_PI); else devf0=1.0/f0; RootSolver<Rayleigh> solver(waveModel); RealValue * point=mode(modeIndex); int n=xCount(); int i,j; for(j=1; j<n && !point[j].isValid(); j++) {} if(j==n) return 0.0; for(i=j+1; i<n && !point[i].isValid(); i++ ) {} if(i==n) return 0.0; double lastSlope=fabs(point[i].value())-fabs(point[j].value()); j=i; for(i++; i<n; i++ ) { if(point[i].isValid()) { double slope=fabs(point[i].value())-fabs(point[j].value()); j=i; if(slope<=0 && lastSlope>0) { double diff; // Found a maximum in the sampled curve that is at distance > 0.01 the ellipticity of the half space if(((diff=(f0-x(i))*devf0)>=0 && diff<misfit) || ((diff=(x(i-2)-f0)*devf0)>=0 && diff<misfit) || (x(i-2)<f0 && f0<x(i))) { // There's a chance to get a better misfit bool badSampling; double xr=refineMax(modeIndex, i, disp, solver, badSampling); //App::stream() << tr("Refined peak frequency=%1").arg(0.5*xr/M_PI) << endl; // The peak is sufficiently refined, adjust the misfit if necessary diff=fabs(xr-f0)*devf0; if(diff<misfit) misfit=diff; if(badSampling) App::stream() << tr("** Warning ** : probably a missing peak, bad initial sampling") << endl; } } lastSlope=slope; } } if(misfit==1e99) { // monotonous curve, evaluate fundamental frequency by V0/4H double roughFreq=2.0 * M_PI * waveModel->model()->roughFundamentalFrequency(); //App::stream() << tr("Rough fundamental frequency=%1").arg(0.5 * roughFreq/M_PI) << endl; if(roughFreq>=x(0) && roughFreq<=x(n-1)) { App::stream() << tr("** Warning ** : estimated rough fundamental frequency inside range") << endl; double midRange=0.5*(x(0)+x(n-1)); if(roughFreq>midRange) misfit=fabs(x(n-1)-f0)*devf0; else misfit=fabs(x(0)-f0)*devf0; } else misfit=fabs(roughFreq-f0)*devf0; } refineSort(); disp.refineSort(); return floor(misfit*1000+0.5)*0.001; // frequency is estimated down to 1e-3, avoid classification within precision limits }
QList< double > QGpCoreWave::Ellipticity::peaks | ( | int | modeIndex, |
Dispersion & | disp, | ||
Rayleigh * | waveModel | ||
) |
Return a list of the frequency of all peaks within the sampled range
References QGpCoreTools::endl(), QGpCoreTools::Value< numberType >::isValid(), QGpCoreWave::ModalStorage::mode(), QGpCoreTools::tr(), QGpCoreTools::Value< numberType >::value(), and QGpCoreWave::ModalStorage::xCount().
Referenced by EllipticityReader::parse().
{ QList<double> peaks; double fac=0.5/M_PI; RootSolver<Rayleigh> solver(waveModel); RealValue * point=mode(modeIndex); int n=xCount(); int i,j; for(j=1; j<n && !point[j].isValid(); j++) {} if(j==n) return peaks; for(i=j+1; i<n && !point[i].isValid(); i++ ) {} if(i==n) return peaks; double lastSlope=point[i].value()-point[j].value(); j=i; for(i++; i<n; i++ ) { if(point[i].isValid()) { double slope=fabs(point[i].value())-fabs(point[j].value()); j=i; if(slope<=0 && lastSlope>0 && lastSlope-slope>1e-10) { // Found a maximum in the sampled curve that is at distance > 0.01 the ellipticity of the half space bool badSampling; double xr=refineMax(modeIndex, i, disp, solver, badSampling); peaks << fac*xr; if(badSampling) App::stream() << tr("** Warning ** : probably a missing peak, bad initial sampling") << endl; } lastSlope=slope; } } return peaks; }