mirror of
https://github.com/OPM/ResInsight.git
synced 2025-02-25 18:55:39 -06:00
Statistics UI is now functional, together with two percentile evaluation methods
Refurbished the calculation algorithm code somewhat p4#: 21134
This commit is contained in:
Jacob Støren
parent
cc8568cca6
commit
75f089d0d4
@ -185,6 +185,12 @@ void RimStatisticsCase::computeStatistics()
|
||||
|
||||
RimStatisticsConfig statisticsConfig;
|
||||
|
||||
statisticsConfig.m_calculatePercentiles = m_calculatePercentiles();
|
||||
statisticsConfig.m_pMaxPos = m_highPercentile();
|
||||
statisticsConfig.m_pMidPos = m_midPercentile();
|
||||
statisticsConfig.m_pMinPos = m_lowPercentile();
|
||||
statisticsConfig.m_pValMethod = m_percentileCalculationType();
|
||||
|
||||
std::vector<size_t> timeStepIndices;
|
||||
for (size_t i = 0; i < timeStepCount; i++)
|
||||
{
|
||||
@ -193,31 +199,53 @@ void RimStatisticsCase::computeStatistics()
|
||||
|
||||
RigCaseData* resultCase = reservoirData();
|
||||
|
||||
QList<QPair<RimDefines::ResultCatType, QString> > resultSpecification;
|
||||
QList<RimStatisticsCaseEvaluator::ResSpec > resultSpecification;
|
||||
|
||||
//resultSpecification.append(qMakePair(RimDefines::DYNAMIC_NATIVE, QString("PRESSURE")));
|
||||
|
||||
|
||||
for(size_t pIdx = 0; pIdx < m_selectedDynamicProperties().size(); ++pIdx)
|
||||
{
|
||||
QStringList resultNames = sourceCases.at(0)->results(RifReaderInterface::MATRIX_RESULTS)->cellResults()->resultNames(RimDefines::DYNAMIC_NATIVE);
|
||||
foreach(QString resultName, resultNames)
|
||||
{
|
||||
resultSpecification.append(qMakePair(RimDefines::DYNAMIC_NATIVE, resultName));
|
||||
}
|
||||
resultSpecification.append(RimStatisticsCaseEvaluator::ResSpec(RifReaderInterface::MATRIX_RESULTS, RimDefines::DYNAMIC_NATIVE, m_selectedDynamicProperties()[pIdx]));
|
||||
}
|
||||
|
||||
for(size_t pIdx = 0; pIdx < m_selectedStaticProperties().size(); ++pIdx)
|
||||
{
|
||||
QStringList resultNames = sourceCases.at(0)->results(RifReaderInterface::MATRIX_RESULTS)->cellResults()->resultNames(RimDefines::STATIC_NATIVE);
|
||||
foreach(QString resultName, resultNames)
|
||||
{
|
||||
resultSpecification.append(qMakePair(RimDefines::STATIC_NATIVE, resultName));
|
||||
}
|
||||
resultSpecification.append(RimStatisticsCaseEvaluator::ResSpec(RifReaderInterface::MATRIX_RESULTS, RimDefines::STATIC_NATIVE, m_selectedStaticProperties()[pIdx]));
|
||||
}
|
||||
|
||||
for(size_t pIdx = 0; pIdx < m_selectedGeneratedProperties().size(); ++pIdx)
|
||||
{
|
||||
resultSpecification.append(RimStatisticsCaseEvaluator::ResSpec(RifReaderInterface::MATRIX_RESULTS, RimDefines::GENERATED, m_selectedGeneratedProperties()[pIdx]));
|
||||
}
|
||||
|
||||
for(size_t pIdx = 0; pIdx < m_selectedInputProperties().size(); ++pIdx)
|
||||
{
|
||||
resultSpecification.append(RimStatisticsCaseEvaluator::ResSpec(RifReaderInterface::MATRIX_RESULTS, RimDefines::INPUT_PROPERTY, m_selectedInputProperties()[pIdx]));
|
||||
}
|
||||
|
||||
for(size_t pIdx = 0; pIdx < m_selectedFractureDynamicProperties().size(); ++pIdx)
|
||||
{
|
||||
resultSpecification.append(RimStatisticsCaseEvaluator::ResSpec(RifReaderInterface::FRACTURE_RESULTS, RimDefines::DYNAMIC_NATIVE, m_selectedFractureDynamicProperties()[pIdx]));
|
||||
}
|
||||
|
||||
for(size_t pIdx = 0; pIdx < m_selectedFractureStaticProperties().size(); ++pIdx)
|
||||
{
|
||||
resultSpecification.append(RimStatisticsCaseEvaluator::ResSpec(RifReaderInterface::FRACTURE_RESULTS, RimDefines::STATIC_NATIVE, m_selectedFractureStaticProperties()[pIdx]));
|
||||
}
|
||||
|
||||
for(size_t pIdx = 0; pIdx < m_selectedFractureGeneratedProperties().size(); ++pIdx)
|
||||
{
|
||||
resultSpecification.append(RimStatisticsCaseEvaluator::ResSpec(RifReaderInterface::FRACTURE_RESULTS, RimDefines::GENERATED, m_selectedFractureGeneratedProperties()[pIdx]));
|
||||
}
|
||||
|
||||
for(size_t pIdx = 0; pIdx < m_selectedFractureInputProperties().size(); ++pIdx)
|
||||
{
|
||||
resultSpecification.append(RimStatisticsCaseEvaluator::ResSpec(RifReaderInterface::FRACTURE_RESULTS, RimDefines::INPUT_PROPERTY, m_selectedFractureInputProperties()[pIdx]));
|
||||
}
|
||||
|
||||
|
||||
RimStatisticsCaseEvaluator stat(sourceCases, timeStepIndices, statisticsConfig, resultCase);
|
||||
stat.evaluateForResults(resultSpecification);
|
||||
|
||||
// Todo: Is this really the time and place to do the following ? JJS
|
||||
|
||||
for (size_t i = 0; i < reservoirViews().size(); i++)
|
||||
{
|
||||
RimReservoirView* reservoirView = reservoirViews()[i];
|
||||
|
||||
@ -30,80 +30,98 @@
|
||||
//--------------------------------------------------------------------------------------------------
|
||||
/// An internal class to do the actual computations
|
||||
//--------------------------------------------------------------------------------------------------
|
||||
class RimStatisticsEvaluator
|
||||
|
||||
void calculateBasicStatistics(const std::vector<double>& values, double* min, double* max, double* range, double* mean, double* dev)
|
||||
{
|
||||
public:
|
||||
RimStatisticsEvaluator(const std::vector<double>& values)
|
||||
: m_values(values),
|
||||
m_min(HUGE_VAL),
|
||||
m_max(-HUGE_VAL),
|
||||
m_mean(HUGE_VAL),
|
||||
m_dev(HUGE_VAL)
|
||||
double m_min(HUGE_VAL);
|
||||
double m_max(-HUGE_VAL);
|
||||
double m_mean(HUGE_VAL);
|
||||
double m_dev(HUGE_VAL);
|
||||
|
||||
double sum = 0.0;
|
||||
double sumSquared = 0.0;
|
||||
|
||||
size_t validValueCount = 0;
|
||||
|
||||
for (size_t i = 0; i < values.size(); i++)
|
||||
{
|
||||
double val = values[i];
|
||||
if (val == HUGE_VAL) continue;
|
||||
|
||||
validValueCount++;
|
||||
|
||||
if (val < m_min) m_min = val;
|
||||
if (val > m_max) m_max = val;
|
||||
|
||||
sum += val;
|
||||
sumSquared += (val * val);
|
||||
}
|
||||
|
||||
|
||||
void getStatistics(double& min, double& max, double& mean, double& dev, double& range)
|
||||
if (validValueCount > 0)
|
||||
{
|
||||
evaluate();
|
||||
m_mean = sum / validValueCount;
|
||||
|
||||
min = m_min;
|
||||
max = m_max;
|
||||
mean = m_mean;
|
||||
dev = m_dev;
|
||||
|
||||
range = m_max - m_min;
|
||||
// http://en.wikipedia.org/wiki/Standard_deviation#Rapid_calculation_methods
|
||||
// Running standard deviation
|
||||
|
||||
double s0 = validValueCount;
|
||||
double s1 = sum;
|
||||
double s2 = sumSquared;
|
||||
|
||||
m_dev = cvf::Math::sqrt( (s0 * s2) - (s1 * s1) ) / s0;
|
||||
}
|
||||
|
||||
private:
|
||||
void evaluate()
|
||||
if (min) *min = m_min;
|
||||
if (max) *max = m_max;
|
||||
if (range) *range = m_max - m_min;
|
||||
|
||||
if (mean) *mean = m_mean;
|
||||
if (dev) *dev = m_dev;
|
||||
}
|
||||
|
||||
|
||||
//--------------------------------------------------------------------------------------------------
|
||||
/// Calculate the percentiles of /a inputValues at the pValPosition percentages using the "Nearest Rank"
|
||||
/// method. This method treats HUGE_VAL as "undefined" values, and ignores these. Will return HUGE_VAL if
|
||||
/// the inputValues does not contain any valid values
|
||||
//--------------------------------------------------------------------------------------------------
|
||||
|
||||
std::vector<double> calculateNearestRankPercentiles(const std::vector<double> & inputValues, const std::vector<double>& pValPositions)
|
||||
{
|
||||
std::vector<double> sortedValues;
|
||||
sortedValues.reserve(inputValues.size());
|
||||
|
||||
for (size_t i = 0; i < inputValues.size(); ++i)
|
||||
{
|
||||
double sum = 0.0;
|
||||
double sumSquared = 0.0;
|
||||
|
||||
size_t validValueCount = 0;
|
||||
|
||||
for (size_t i = 0; i < m_values.size(); i++)
|
||||
if (inputValues[i] != HUGE_VAL)
|
||||
{
|
||||
double val = m_values[i];
|
||||
if (val == HUGE_VAL) continue;
|
||||
|
||||
validValueCount++;
|
||||
|
||||
if (val < m_min) m_min = val;
|
||||
if (val > m_max) m_max = val;
|
||||
|
||||
sum += val;
|
||||
sumSquared += (val * val);
|
||||
}
|
||||
|
||||
if (validValueCount > 0)
|
||||
{
|
||||
m_mean = sum / validValueCount;
|
||||
|
||||
|
||||
// http://en.wikipedia.org/wiki/Standard_deviation#Rapid_calculation_methods
|
||||
// Running standard deviation
|
||||
|
||||
double s0 = validValueCount;
|
||||
double s1 = sum;
|
||||
double s2 = sumSquared;
|
||||
|
||||
m_dev = cvf::Math::sqrt( (s0 * s2) - (s1 * s1) ) / s0;
|
||||
sortedValues.push_back(inputValues[i]);
|
||||
}
|
||||
}
|
||||
|
||||
std::sort(sortedValues.begin(), sortedValues.end());
|
||||
|
||||
private:
|
||||
const std::vector<double>& m_values;
|
||||
std::vector<double> percentiles(pValPositions.size(), HUGE_VAL);
|
||||
if (sortedValues.size())
|
||||
{
|
||||
for (size_t i = 0; i < pValPositions.size(); ++i)
|
||||
{
|
||||
double pVal = HUGE_VAL;
|
||||
|
||||
double m_min;
|
||||
double m_max;
|
||||
double m_mean;
|
||||
double m_dev;
|
||||
size_t pValIndex = static_cast<size_t>(sortedValues.size() * abs(pValPositions[i]) / 100);
|
||||
pValIndex += 1;
|
||||
|
||||
if (pValIndex >= sortedValues.size() ) pValIndex = sortedValues.size() - 1;
|
||||
|
||||
pVal = sortedValues[pValIndex];
|
||||
percentiles[i] = pVal;
|
||||
}
|
||||
}
|
||||
|
||||
return percentiles;
|
||||
};
|
||||
|
||||
|
||||
//--------------------------------------------------------------------------------------------------
|
||||
///
|
||||
//--------------------------------------------------------------------------------------------------
|
||||
@ -135,29 +153,28 @@ QString createResultNameMax(const QString& resultName) { return resultName + "_
|
||||
QString createResultNameMean(const QString& resultName) { return resultName + "_MEAN"; }
|
||||
QString createResultNameDev(const QString& resultName) { return resultName + "_DEV"; }
|
||||
QString createResultNameRange(const QString& resultName) { return resultName + "_RANGE"; }
|
||||
QString createResultNamePVal(const QString& resultName, double pValPos) { return resultName + "_P_" + QString::number(pValPos); }
|
||||
|
||||
|
||||
//--------------------------------------------------------------------------------------------------
|
||||
///
|
||||
//--------------------------------------------------------------------------------------------------
|
||||
void RimStatisticsCaseEvaluator::buildSourceMetaData(RimDefines::ResultCatType resultType, const QString& resultName)
|
||||
void RimStatisticsCaseEvaluator::buildSourceMetaData(RifReaderInterface::PorosityModelResultType poroModel, RimDefines::ResultCatType resultType, const QString& resultName)
|
||||
{
|
||||
if (m_sourceCases.size() == 0) return;
|
||||
|
||||
std::vector<QDateTime> timeStepDates = m_sourceCases[0]->results(RifReaderInterface::MATRIX_RESULTS)->cellResults()->timeStepDates(0);
|
||||
std::vector<QDateTime> timeStepDates = m_sourceCases[0]->results(poroModel)->cellResults()->timeStepDates(0);
|
||||
|
||||
for (size_t caseIdx = 1; caseIdx < m_sourceCases.size(); caseIdx++)
|
||||
{
|
||||
RigCaseData* eclipseCase = m_sourceCases.at(caseIdx)->reservoirData();
|
||||
|
||||
RimReservoirCellResultsStorage* matrixResults = m_sourceCases[caseIdx]->results(RifReaderInterface::MATRIX_RESULTS);
|
||||
size_t scalarResultIndex = matrixResults->findOrLoadScalarResult(resultType, resultName);
|
||||
RimReservoirCellResultsStorage* cellResultsStorage = m_sourceCases[caseIdx]->results(poroModel);
|
||||
size_t scalarResultIndex = cellResultsStorage->findOrLoadScalarResult(resultType, resultName);
|
||||
if (scalarResultIndex == cvf::UNDEFINED_SIZE_T)
|
||||
{
|
||||
size_t scalarResultIndex = matrixResults->cellResults()->addEmptyScalarResult(resultType, resultName, false);
|
||||
matrixResults->cellResults()->setTimeStepDates(scalarResultIndex, timeStepDates);
|
||||
size_t scalarResultIndex = cellResultsStorage->cellResults()->addEmptyScalarResult(resultType, resultName, false);
|
||||
cellResultsStorage->cellResults()->setTimeStepDates(scalarResultIndex, timeStepDates);
|
||||
|
||||
std::vector< std::vector<double> >& dataValues = matrixResults->cellResults()->cellScalarResults(scalarResultIndex);
|
||||
std::vector< std::vector<double> >& dataValues = cellResultsStorage->cellResults()->cellScalarResults(scalarResultIndex);
|
||||
dataValues.resize(timeStepDates.size());
|
||||
}
|
||||
}
|
||||
@ -166,222 +183,234 @@ void RimStatisticsCaseEvaluator::buildSourceMetaData(RimDefines::ResultCatType r
|
||||
//--------------------------------------------------------------------------------------------------
|
||||
///
|
||||
//--------------------------------------------------------------------------------------------------
|
||||
void RimStatisticsCaseEvaluator::evaluateForResults(const QList<QPair<RimDefines::ResultCatType, QString> >& resultSpecification)
|
||||
void RimStatisticsCaseEvaluator::evaluateForResults(const QList<ResSpec>& resultSpecification)
|
||||
{
|
||||
CVF_ASSERT(m_destinationCase);
|
||||
|
||||
size_t activeMatrixCellCount = m_destinationCase->activeCellInfo(RifReaderInterface::MATRIX_RESULTS)->globalActiveCellCount();
|
||||
RigCaseCellResultsData* matrixResults = m_destinationCase->results(RifReaderInterface::MATRIX_RESULTS);
|
||||
|
||||
// First build the destination result data structures to receive the statistics
|
||||
|
||||
for (int i = 0; i < resultSpecification.size(); i++)
|
||||
{
|
||||
RimDefines::ResultCatType resultType = resultSpecification[i].first;
|
||||
QString resultName = resultSpecification[i].second;
|
||||
RifReaderInterface::PorosityModelResultType poroModel = resultSpecification[i].m_poroModel;
|
||||
RimDefines::ResultCatType resultType = resultSpecification[i].m_resType;
|
||||
QString resultName = resultSpecification[i].m_resVarName;
|
||||
|
||||
size_t activeCellCount = m_destinationCase->activeCellInfo(poroModel)->globalActiveCellCount();
|
||||
RigCaseCellResultsData* destCellResultsData = m_destinationCase->results(poroModel);
|
||||
|
||||
// Special handling if SOIL is asked for
|
||||
// Build SGAS/SWAT meta data, SOIL is automatically generated as part of RigCaseCellResultsData::findOrLoadScalarResultForTimeStep
|
||||
if (resultName.toUpper() == "SOIL")
|
||||
{
|
||||
size_t swatIndex = m_sourceCases.at(0)->results(RifReaderInterface::MATRIX_RESULTS)->cellResults()->findScalarResultIndex(resultType, "SWAT");
|
||||
size_t swatIndex = m_sourceCases.at(0)->results(poroModel)->cellResults()->findScalarResultIndex(resultType, "SWAT");
|
||||
if (swatIndex != cvf::UNDEFINED_SIZE_T)
|
||||
{
|
||||
buildSourceMetaData(resultType, "SWAT");
|
||||
buildSourceMetaData(poroModel, resultType, "SWAT");
|
||||
}
|
||||
|
||||
size_t sgasIndex = m_sourceCases.at(0)->results(RifReaderInterface::MATRIX_RESULTS)->cellResults()->findScalarResultIndex(resultType, "SGAS");
|
||||
size_t sgasIndex = m_sourceCases.at(0)->results(poroModel)->cellResults()->findScalarResultIndex(resultType, "SGAS");
|
||||
if (sgasIndex != cvf::UNDEFINED_SIZE_T)
|
||||
{
|
||||
buildSourceMetaData(resultType, "SGAS");
|
||||
buildSourceMetaData(poroModel, resultType, "SGAS");
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
// Meta info is loaded from disk for first case only
|
||||
// Build metadata for all other source cases
|
||||
buildSourceMetaData(resultType, resultName);
|
||||
// Meta info is loaded from disk for first case only
|
||||
// Build metadata for all other source cases
|
||||
buildSourceMetaData(poroModel, resultType, resultName);
|
||||
}
|
||||
|
||||
QString minResultName = createResultNameMin(resultName);
|
||||
QString maxResultName = createResultNameMax(resultName);
|
||||
QString meanResultName = createResultNameMean(resultName);
|
||||
QString devResultName = createResultNameDev(resultName);
|
||||
QString rangeResultName = createResultNameRange(resultName);
|
||||
// Create new result data structures to contain the statistical values
|
||||
std::vector<QString> statisticalResultNames;
|
||||
|
||||
if (activeMatrixCellCount > 0)
|
||||
statisticalResultNames.push_back(createResultNameMin(resultName));
|
||||
statisticalResultNames.push_back(createResultNameMax(resultName));
|
||||
statisticalResultNames.push_back(createResultNameMean(resultName));
|
||||
statisticalResultNames.push_back(createResultNameDev(resultName));
|
||||
statisticalResultNames.push_back(createResultNameRange(resultName));
|
||||
|
||||
if (m_statisticsConfig.m_calculatePercentiles)
|
||||
{
|
||||
addNamedResult(matrixResults, resultType, minResultName, activeMatrixCellCount);
|
||||
addNamedResult(matrixResults, resultType, maxResultName, activeMatrixCellCount);
|
||||
addNamedResult(matrixResults, resultType, meanResultName, activeMatrixCellCount);
|
||||
addNamedResult(matrixResults, resultType, devResultName, activeMatrixCellCount);
|
||||
addNamedResult(matrixResults, resultType, rangeResultName, activeMatrixCellCount);
|
||||
statisticalResultNames.push_back(createResultNamePVal(resultName, m_statisticsConfig.m_pMinPos));
|
||||
statisticalResultNames.push_back(createResultNamePVal(resultName, m_statisticsConfig.m_pMidPos));
|
||||
statisticalResultNames.push_back(createResultNamePVal(resultName, m_statisticsConfig.m_pMaxPos));
|
||||
}
|
||||
|
||||
if (activeCellCount > 0)
|
||||
{
|
||||
for (size_t i = 0; i < statisticalResultNames.size(); ++i)
|
||||
{
|
||||
addNamedResult(destCellResultsData, resultType, statisticalResultNames[i], activeCellCount);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
if (activeMatrixCellCount > 0)
|
||||
// Start the loop that calculates the statistics
|
||||
|
||||
caf::ProgressInfo progressInfo(m_timeStepIndices.size(), "Computing Statistics");
|
||||
|
||||
for (size_t timeIndicesIdx = 0; timeIndicesIdx < m_timeStepIndices.size(); timeIndicesIdx++)
|
||||
{
|
||||
caf::ProgressInfo info(m_timeStepIndices.size(), "Computing Statistics");
|
||||
size_t timeStepIdx = m_timeStepIndices[timeIndicesIdx];
|
||||
|
||||
for (size_t timeIndicesIdx = 0; timeIndicesIdx < m_timeStepIndices.size(); timeIndicesIdx++)
|
||||
for (size_t gridIdx = 0; gridIdx < m_destinationCase->gridCount(); gridIdx++)
|
||||
{
|
||||
size_t timeStepIdx = m_timeStepIndices[timeIndicesIdx];
|
||||
RigGridBase* grid = m_destinationCase->grid(gridIdx);
|
||||
|
||||
size_t gridCount = 0;
|
||||
for (size_t gridIdx = 0; gridIdx < m_destinationCase->gridCount(); gridIdx++)
|
||||
for (int i = 0; i < resultSpecification.size(); i++)
|
||||
{
|
||||
RigGridBase* grid = m_destinationCase->grid(gridIdx);
|
||||
RifReaderInterface::PorosityModelResultType poroModel = resultSpecification[i].m_poroModel;
|
||||
RimDefines::ResultCatType resultType = resultSpecification[i].m_resType;
|
||||
QString resultName = resultSpecification[i].m_resVarName;
|
||||
|
||||
for (int i = 0; i < resultSpecification.size(); i++)
|
||||
size_t activeCellCount = m_destinationCase->activeCellInfo(poroModel)->globalActiveCellCount();
|
||||
|
||||
if (activeCellCount == 0) continue;
|
||||
|
||||
RigCaseCellResultsData* destCellResultsData = m_destinationCase->results(poroModel);
|
||||
|
||||
size_t dataAccessTimeStepIndex = timeStepIdx;
|
||||
|
||||
// Always evaluate statistics once, and always use time step index zero
|
||||
if (resultType == RimDefines::STATIC_NATIVE)
|
||||
{
|
||||
RimDefines::ResultCatType resultType = resultSpecification[i].first;
|
||||
QString resultName = resultSpecification[i].second;
|
||||
if (timeIndicesIdx > 0) continue;
|
||||
|
||||
size_t dataAccessTimeStepIndex = timeStepIdx;
|
||||
dataAccessTimeStepIndex = 0;
|
||||
}
|
||||
|
||||
// Always evaluate statistics once, and always use time step index zero
|
||||
if (resultType == RimDefines::STATIC_NATIVE)
|
||||
// Build data access objects for source scalar results
|
||||
|
||||
cvf::Collection<cvf::StructGridScalarDataAccess> sourceDataAccessList;
|
||||
for (size_t caseIdx = 0; caseIdx < m_sourceCases.size(); caseIdx++)
|
||||
{
|
||||
RimCase* eclipseCase = m_sourceCases.at(caseIdx);
|
||||
|
||||
size_t scalarResultIndex = eclipseCase->results(poroModel)->findOrLoadScalarResultForTimeStep(resultType, resultName, dataAccessTimeStepIndex);
|
||||
|
||||
cvf::ref<cvf::StructGridScalarDataAccess> dataAccessObject = eclipseCase->reservoirData()->dataAccessObject(grid, poroModel, dataAccessTimeStepIndex, scalarResultIndex);
|
||||
if (dataAccessObject.notNull())
|
||||
{
|
||||
if (timeIndicesIdx > 0) continue;
|
||||
|
||||
dataAccessTimeStepIndex = 0;
|
||||
sourceDataAccessList.push_back(dataAccessObject.p());
|
||||
}
|
||||
}
|
||||
|
||||
// Build data access objects for source scalar results
|
||||
cvf::Collection<cvf::StructGridScalarDataAccess> dataAccesObjectList;
|
||||
for (size_t caseIdx = 0; caseIdx < m_sourceCases.size(); caseIdx++)
|
||||
// Build data access objects for destination scalar results
|
||||
// Find the created result container, if any, and put its dataAccessObject into the enum indexed destination collection
|
||||
|
||||
cvf::Collection<cvf::StructGridScalarDataAccess> destinationDataAccessList;
|
||||
std::vector<QString> statisticalResultNames(STAT_PARAM_COUNT);
|
||||
|
||||
statisticalResultNames[MIN] = createResultNameMin(resultName);
|
||||
statisticalResultNames[MAX] = createResultNameMax(resultName);
|
||||
statisticalResultNames[RANGE] = createResultNameRange(resultName);
|
||||
statisticalResultNames[MEAN] = createResultNameMean(resultName);
|
||||
statisticalResultNames[STDEV] = createResultNameDev(resultName);
|
||||
statisticalResultNames[PMIN] = createResultNamePVal(resultName, m_statisticsConfig.m_pMinPos);
|
||||
statisticalResultNames[PMID] = createResultNamePVal(resultName, m_statisticsConfig.m_pMidPos);
|
||||
statisticalResultNames[PMAX] = createResultNamePVal(resultName, m_statisticsConfig.m_pMaxPos);
|
||||
|
||||
for (size_t stIdx = 0; stIdx < statisticalResultNames.size(); ++stIdx)
|
||||
{
|
||||
size_t scalarResultIndex = destCellResultsData->findScalarResultIndex(resultType, statisticalResultNames[stIdx]);
|
||||
if (scalarResultIndex != cvf::UNDEFINED_SIZE_T)
|
||||
{
|
||||
RimCase* eclipseCase = m_sourceCases.at(caseIdx);
|
||||
|
||||
size_t scalarResultIndex = eclipseCase->results(RifReaderInterface::MATRIX_RESULTS)->findOrLoadScalarResultForTimeStep(resultType, resultName, dataAccessTimeStepIndex);
|
||||
|
||||
cvf::ref<cvf::StructGridScalarDataAccess> dataAccessObject = eclipseCase->reservoirData()->dataAccessObject(grid, RifReaderInterface::MATRIX_RESULTS, dataAccessTimeStepIndex, scalarResultIndex);
|
||||
if (dataAccessObject.notNull())
|
||||
{
|
||||
dataAccesObjectList.push_back(dataAccessObject.p());
|
||||
}
|
||||
destinationDataAccessList.push_back(m_destinationCase->dataAccessObject(grid, poroModel, dataAccessTimeStepIndex, scalarResultIndex).p());
|
||||
}
|
||||
|
||||
|
||||
// Build data access objects form destination scalar results
|
||||
cvf::ref<cvf::StructGridScalarDataAccess> dataAccessObjectMin = NULL;
|
||||
cvf::ref<cvf::StructGridScalarDataAccess> dataAccessObjectMax = NULL;
|
||||
cvf::ref<cvf::StructGridScalarDataAccess> dataAccessObjectMean = NULL;
|
||||
cvf::ref<cvf::StructGridScalarDataAccess> dataAccessObjectDev = NULL;
|
||||
cvf::ref<cvf::StructGridScalarDataAccess> dataAccessObjectRange = NULL;
|
||||
|
||||
else
|
||||
{
|
||||
size_t scalarResultIndex = matrixResults->findScalarResultIndex(resultType, createResultNameMin(resultName));
|
||||
if (scalarResultIndex != cvf::UNDEFINED_SIZE_T)
|
||||
{
|
||||
dataAccessObjectMin = m_destinationCase->dataAccessObject(grid, RifReaderInterface::MATRIX_RESULTS, dataAccessTimeStepIndex, scalarResultIndex);
|
||||
}
|
||||
destinationDataAccessList.push_back(NULL);
|
||||
}
|
||||
}
|
||||
|
||||
// Loop over the cells in the grid, get the case values, and calculate the cell statistics
|
||||
|
||||
for (size_t cellIdx = 0; cellIdx < grid->cellCount(); cellIdx++)
|
||||
{
|
||||
|
||||
size_t globalGridCellIdx = grid->globalGridCellIndex(cellIdx);
|
||||
if (m_destinationCase->activeCellInfo(poroModel)->isActive(globalGridCellIdx))
|
||||
{
|
||||
size_t scalarResultIndex = matrixResults->findScalarResultIndex(resultType, createResultNameMax(resultName));
|
||||
if (scalarResultIndex != cvf::UNDEFINED_SIZE_T)
|
||||
{
|
||||
dataAccessObjectMax = m_destinationCase->dataAccessObject(grid, RifReaderInterface::MATRIX_RESULTS, dataAccessTimeStepIndex, scalarResultIndex);
|
||||
}
|
||||
}
|
||||
// Extract the cell values from each of the cases and assemble them into one vector
|
||||
|
||||
{
|
||||
size_t scalarResultIndex = matrixResults->findScalarResultIndex(resultType, createResultNameMean(resultName));
|
||||
if (scalarResultIndex != cvf::UNDEFINED_SIZE_T)
|
||||
{
|
||||
dataAccessObjectMean = m_destinationCase->dataAccessObject(grid, RifReaderInterface::MATRIX_RESULTS, dataAccessTimeStepIndex, scalarResultIndex);
|
||||
}
|
||||
}
|
||||
std::vector<double> values(sourceDataAccessList.size(), HUGE_VAL);
|
||||
|
||||
{
|
||||
size_t scalarResultIndex = matrixResults->findScalarResultIndex(resultType, createResultNameDev(resultName));
|
||||
if (scalarResultIndex != cvf::UNDEFINED_SIZE_T)
|
||||
bool foundAnyValidValues = false;
|
||||
for (size_t caseIdx = 0; caseIdx < sourceDataAccessList.size(); caseIdx++)
|
||||
{
|
||||
dataAccessObjectDev = m_destinationCase->dataAccessObject(grid, RifReaderInterface::MATRIX_RESULTS, dataAccessTimeStepIndex, scalarResultIndex);
|
||||
}
|
||||
}
|
||||
double val = sourceDataAccessList.at(caseIdx)->cellScalar(cellIdx);
|
||||
values[caseIdx] = val;
|
||||
|
||||
{
|
||||
size_t scalarResultIndex = matrixResults->findScalarResultIndex(resultType, createResultNameRange(resultName));
|
||||
if (scalarResultIndex != cvf::UNDEFINED_SIZE_T)
|
||||
{
|
||||
dataAccessObjectRange = m_destinationCase->dataAccessObject(grid, RifReaderInterface::MATRIX_RESULTS, dataAccessTimeStepIndex, scalarResultIndex);
|
||||
}
|
||||
}
|
||||
|
||||
double min, max, mean, dev, range;
|
||||
for (size_t cellIdx = 0; cellIdx < grid->cellCount(); cellIdx++)
|
||||
{
|
||||
std::vector<double> values(dataAccesObjectList.size(), HUGE_VAL);
|
||||
|
||||
size_t globalGridCellIdx = grid->globalGridCellIndex(cellIdx);
|
||||
if (m_destinationCase->activeCellInfo(RifReaderInterface::MATRIX_RESULTS)->isActive(globalGridCellIdx))
|
||||
{
|
||||
bool foundAnyValidValues = false;
|
||||
for (size_t caseIdx = 0; caseIdx < dataAccesObjectList.size(); caseIdx++)
|
||||
if (val != HUGE_VAL)
|
||||
{
|
||||
double val = dataAccesObjectList.at(caseIdx)->cellScalar(cellIdx);
|
||||
values[caseIdx] = val;
|
||||
foundAnyValidValues = true;
|
||||
}
|
||||
}
|
||||
|
||||
if (val != HUGE_VAL)
|
||||
// Do the real statistics calculations
|
||||
std::vector<double> statParams(STAT_PARAM_COUNT, HUGE_VAL);
|
||||
|
||||
if (foundAnyValidValues)
|
||||
{
|
||||
calculateBasicStatistics(values, &statParams[MIN], &statParams[MAX], &statParams[RANGE], &statParams[MEAN], &statParams[STDEV]);
|
||||
|
||||
// Calculate percentiles
|
||||
if (m_statisticsConfig.m_calculatePercentiles )
|
||||
{
|
||||
if (m_statisticsConfig.m_pValMethod == RimStatisticsCase::EXACT)
|
||||
{
|
||||
foundAnyValidValues = true;
|
||||
std::vector<double> pValPoss;
|
||||
pValPoss.push_back(m_statisticsConfig.m_pMinPos);
|
||||
pValPoss.push_back(m_statisticsConfig.m_pMidPos);
|
||||
pValPoss.push_back(m_statisticsConfig.m_pMaxPos);
|
||||
std::vector<double> pVals = calculateNearestRankPercentiles(values, pValPoss);
|
||||
statParams[PMIN] = pVals[0];
|
||||
statParams[PMID] = pVals[1];
|
||||
statParams[PMAX] = pVals[2];
|
||||
}
|
||||
else if (m_statisticsConfig.m_pValMethod == RimStatisticsCase::HISTOGRAM_ESTIMATED)
|
||||
{
|
||||
std::vector<size_t> histogram;
|
||||
RigHistogramCalculator histCalc(statParams[MIN], statParams[MAX], 100, &histogram);
|
||||
histCalc.addData(values);
|
||||
statParams[PMIN] = histCalc.calculatePercentil(m_statisticsConfig.m_pMinPos);
|
||||
statParams[PMID] = histCalc.calculatePercentil(m_statisticsConfig.m_pMidPos);
|
||||
statParams[PMAX] = histCalc.calculatePercentil(m_statisticsConfig.m_pMaxPos);
|
||||
}
|
||||
else
|
||||
{
|
||||
CVF_ASSERT(false);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
min = HUGE_VAL;
|
||||
max = HUGE_VAL;
|
||||
mean = HUGE_VAL;
|
||||
dev = HUGE_VAL;
|
||||
range = HUGE_VAL;
|
||||
|
||||
if (foundAnyValidValues)
|
||||
{
|
||||
RimStatisticsEvaluator stat(values);
|
||||
stat.getStatistics(min, max, mean, dev, range);
|
||||
}
|
||||
// Set the results into the results data structures
|
||||
|
||||
if (dataAccessObjectMin.notNull())
|
||||
for (size_t stIdx = 0; stIdx < statParams.size(); ++stIdx)
|
||||
{
|
||||
if (destinationDataAccessList[stIdx].notNull())
|
||||
{
|
||||
dataAccessObjectMin->setCellScalar(cellIdx, min);
|
||||
}
|
||||
|
||||
if (dataAccessObjectMax.notNull())
|
||||
{
|
||||
dataAccessObjectMax->setCellScalar(cellIdx, max);
|
||||
}
|
||||
|
||||
if (dataAccessObjectMean.notNull())
|
||||
{
|
||||
dataAccessObjectMean->setCellScalar(cellIdx, mean);
|
||||
}
|
||||
|
||||
if (dataAccessObjectDev.notNull())
|
||||
{
|
||||
dataAccessObjectDev->setCellScalar(cellIdx, dev);
|
||||
}
|
||||
|
||||
if (dataAccessObjectRange.notNull())
|
||||
{
|
||||
dataAccessObjectRange->setCellScalar(cellIdx, range);
|
||||
destinationDataAccessList[stIdx]->setCellScalar(cellIdx, statParams[stIdx]);
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
for (size_t caseIdx = 0; caseIdx < m_sourceCases.size(); caseIdx++)
|
||||
{
|
||||
RimCase* eclipseCase = m_sourceCases.at(caseIdx);
|
||||
|
||||
// When one time step is completed, close all result files.
|
||||
// Microsoft note: On Windows, the maximum number of files open at the same time is 512
|
||||
// http://msdn.microsoft.com/en-us/library/kdfaxaay%28vs.71%29.aspx
|
||||
//
|
||||
eclipseCase->results(RifReaderInterface::MATRIX_RESULTS)->readerInterface()->close();
|
||||
}
|
||||
|
||||
info.setProgress(timeIndicesIdx);
|
||||
}
|
||||
|
||||
// When one time step is completed, close all result files.
|
||||
// Microsoft note: On Windows, the maximum number of files open at the same time is 512
|
||||
// http://msdn.microsoft.com/en-us/library/kdfaxaay%28vs.71%29.aspx
|
||||
|
||||
for (size_t caseIdx = 0; caseIdx < m_sourceCases.size(); caseIdx++)
|
||||
{
|
||||
RimCase* eclipseCase = m_sourceCases.at(caseIdx);
|
||||
eclipseCase->results(RifReaderInterface::MATRIX_RESULTS)->readerInterface()->close();
|
||||
eclipseCase->results(RifReaderInterface::FRACTURE_RESULTS)->readerInterface()->close();
|
||||
}
|
||||
|
||||
progressInfo.setProgress(timeIndicesIdx);
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
@ -26,6 +26,7 @@
|
||||
|
||||
#include <QPair>
|
||||
#include "RimDefines.h"
|
||||
#include "RimStatisticsCase.h"
|
||||
|
||||
class RimCase;
|
||||
class RigCaseData;
|
||||
@ -36,18 +37,20 @@ class RimStatisticsConfig
|
||||
{
|
||||
public:
|
||||
RimStatisticsConfig()
|
||||
: m_min(true),
|
||||
m_max(true),
|
||||
m_mean(true),
|
||||
m_stdDev(true)
|
||||
: m_calculatePercentiles(true),
|
||||
m_pMinPos(10.0),
|
||||
m_pMidPos(50.0),
|
||||
m_pMaxPos(90.0),
|
||||
m_pValMethod(RimStatisticsCase::EXACT)
|
||||
{
|
||||
}
|
||||
|
||||
public:
|
||||
bool m_min;
|
||||
bool m_max;
|
||||
bool m_mean;
|
||||
bool m_stdDev;
|
||||
bool m_calculatePercentiles;
|
||||
double m_pMinPos;
|
||||
double m_pMidPos;
|
||||
double m_pMaxPos;
|
||||
RimStatisticsCase::PercentileCalcType m_pValMethod;
|
||||
};
|
||||
|
||||
|
||||
@ -59,21 +62,34 @@ public:
|
||||
const RimStatisticsConfig& statisticsConfig,
|
||||
RigCaseData* destinationCase);
|
||||
|
||||
struct ResSpec
|
||||
{
|
||||
ResSpec() : m_resType(RimDefines::DYNAMIC_NATIVE), m_poroModel(RifReaderInterface::MATRIX_RESULTS) {}
|
||||
ResSpec( RifReaderInterface::PorosityModelResultType poroModel,
|
||||
RimDefines::ResultCatType resType,
|
||||
QString resVarName) : m_poroModel(poroModel), m_resType(resType), m_resVarName(resVarName) {}
|
||||
|
||||
void evaluateForResults(const QList<QPair<RimDefines::ResultCatType, QString> >& resultSpecification);
|
||||
RifReaderInterface::PorosityModelResultType m_poroModel;
|
||||
RimDefines::ResultCatType m_resType;
|
||||
QString m_resVarName;
|
||||
};
|
||||
|
||||
void evaluateForResults(const QList<ResSpec >& resultSpecification);
|
||||
|
||||
void debugOutput(RimDefines::ResultCatType resultType, const QString& resultName, size_t timeStepIdx);
|
||||
|
||||
private:
|
||||
void addNamedResult(RigCaseCellResultsData* cellResults, RimDefines::ResultCatType resultType, const QString& resultName, size_t activeCellCount);
|
||||
void buildSourceMetaData(RimDefines::ResultCatType resultType, const QString& resultName);
|
||||
void buildSourceMetaData(RifReaderInterface::PorosityModelResultType poroModel, RimDefines::ResultCatType resultType, const QString& resultName);
|
||||
|
||||
enum StatisticsParamType { MIN, MAX, RANGE, MEAN, STDEV, PMIN, PMID, PMAX, STAT_PARAM_COUNT };
|
||||
|
||||
private:
|
||||
std::vector<RimCase*> m_sourceCases;
|
||||
std::vector<size_t> m_timeStepIndices;
|
||||
std::vector<size_t> m_timeStepIndices;
|
||||
|
||||
size_t m_globalCellCount;
|
||||
RimStatisticsConfig m_statisticsConfig;
|
||||
RigCaseData* m_destinationCase;
|
||||
size_t m_globalCellCount;
|
||||
RimStatisticsConfig m_statisticsConfig;
|
||||
RigCaseData* m_destinationCase;
|
||||
};
|
||||
|
||||
|
||||
@ -180,7 +180,9 @@ const std::vector<size_t>& RigCaseCellResultsData::cellScalarValuesHistogram(siz
|
||||
//--------------------------------------------------------------------------------------------------
|
||||
void RigCaseCellResultsData::p10p90CellScalarValues(size_t scalarResultIndex, double& p10, double& p90)
|
||||
{
|
||||
const std::vector<size_t>& histogr = cellScalarValuesHistogram( scalarResultIndex);
|
||||
// First make sure they are calculated
|
||||
const std::vector<size_t>& histogr = cellScalarValuesHistogram( scalarResultIndex);
|
||||
// Then return them
|
||||
p10 = m_p10p90[scalarResultIndex].first;
|
||||
p90 = m_p10p90[scalarResultIndex].second;
|
||||
}
|
||||
|
||||
Loading…
Reference in New Issue
Block a user