(#606) (#607) Calculate visible cells statistics for geomech

Refactored some in 3D info regarding geomech cases.
This commit is contained in:
Jacob Støren
2015-11-06 10:18:55 +01:00
parent 2adb9279ce
commit 0bd4f4a8f9
11 changed files with 509 additions and 174 deletions

View File

@@ -27,6 +27,8 @@ add_library( ${PROJECT_NAME}
RigFemScalarResultFrames.cpp
RigFemNativeStatCalc.h
RigFemNativeStatCalc.cpp
RigFemNativeVisibleCellsStatCalc.h
RigFemNativeVisibleCellsStatCalc.cpp
RigFemFaceComparator.h
RigFemPartGrid.h
RigFemPartGrid.cpp

View File

@@ -38,7 +38,7 @@ RigFemNativeStatCalc::RigFemNativeStatCalc(RigFemPartResultsCollection* femResul
//--------------------------------------------------------------------------------------------------
void RigFemNativeStatCalc::minMaxCellScalarValues(size_t timeStepIndex, double& min, double& max)
{
for (int pIdx = 0; pIdx < static_cast<int>(m_resultsData->m_femPartResults.size()); ++pIdx)
for (int pIdx = 0; pIdx < m_resultsData->partCount(); ++pIdx)
{
const std::vector<float>& values = m_resultsData->resultValues(m_resVarAddr, pIdx, (int)timeStepIndex);
@@ -68,7 +68,7 @@ void RigFemNativeStatCalc::minMaxCellScalarValues(size_t timeStepIndex, double&
//--------------------------------------------------------------------------------------------------
void RigFemNativeStatCalc::posNegClosestToZero(size_t timeStepIndex, double& pos, double& neg)
{
for (int pIdx = 0; pIdx < static_cast<int>(m_resultsData->m_femPartResults.size()); ++pIdx)
for (int pIdx = 0; pIdx < m_resultsData->partCount(); ++pIdx)
{
const std::vector<float>& values = m_resultsData->resultValues(m_resVarAddr, pIdx, (int)timeStepIndex);
@@ -99,7 +99,7 @@ void RigFemNativeStatCalc::posNegClosestToZero(size_t timeStepIndex, double& pos
void RigFemNativeStatCalc::valueSumAndSampleCount(size_t timeStepIndex, double& valueSum, size_t& sampleCount)
{
int tsIdx = static_cast<int>(timeStepIndex);
int partCount = static_cast<int>(m_resultsData->m_femPartResults.size());
int partCount = m_resultsData->partCount();
for (int pIdx = 0; pIdx < partCount; ++pIdx)
{
@@ -128,7 +128,7 @@ void RigFemNativeStatCalc::valueSumAndSampleCount(size_t timeStepIndex, double&
//--------------------------------------------------------------------------------------------------
void RigFemNativeStatCalc::addDataToHistogramCalculator(size_t timeStepIndex, RigHistogramCalculator& histogramCalculator)
{
int partCount = static_cast<int>(m_resultsData->m_femPartResults.size());
int partCount = m_resultsData->partCount();
for (int pIdx = 0; pIdx < partCount; ++pIdx)
{
const std::vector<float>& values = m_resultsData->resultValues(m_resVarAddr, pIdx, static_cast<int>(timeStepIndex));

View File

@@ -0,0 +1,166 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2015- Statoil ASA
// Copyright (C) 2015- Ceetron Solutions AS
//
// ResInsight 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.
//
// ResInsight 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 at <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#include "RigFemNativeVisibleCellsStatCalc.h"
#include "RigFemScalarResultFrames.h"
#include "RigFemPartResultsCollection.h"
#include <math.h>
#include "RigStatisticsMath.h"
#include "RigGeoMechCaseData.h"
#include "RigFemPartCollection.h"
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigFemNativeVisibleCellsStatCalc::RigFemNativeVisibleCellsStatCalc(RigGeoMechCaseData* femCase,
const RigFemResultAddress& resVarAddr,
const cvf::UByteArray* cellVisibilities)
: m_caseData(femCase), m_resVarAddr(resVarAddr), m_cellVisibilities(cellVisibilities)
{
m_resultsData = femCase->femPartResults();
}
class MinMaxAccumulator
{
public:
MinMaxAccumulator(double initMin, double initMax): max(initMax), min(initMin) {}
void addValue(double value)
{
if (value == HUGE_VAL) // TODO
{
return;
}
if (value < min)
{
min = value;
}
if (value > max)
{
max = value;
}
}
double max;
double min;
};
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigFemNativeVisibleCellsStatCalc::minMaxCellScalarValues(size_t timeStepIndex, double& min, double& max)
{
MinMaxAccumulator acc(min, max);
traverseElementNodes(acc, timeStepIndex);
min = acc.min;
max = acc.max;
}
class PosNegAccumulator
{
public:
PosNegAccumulator(double initPos, double initNeg): pos(initPos), neg(initNeg) {}
void addValue(double value)
{
if (value == HUGE_VAL)
{
return;
}
if (value < pos && value > 0)
{
pos = value;
}
if (value > neg && value < 0)
{
neg = value;
}
}
double pos;
double neg;
};
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigFemNativeVisibleCellsStatCalc::posNegClosestToZero(size_t timeStepIndex, double& pos, double& neg)
{
PosNegAccumulator acc(pos, neg);
traverseElementNodes(acc, timeStepIndex);
pos = acc.pos;
neg = acc.neg;
}
class SumCountAccumulator
{
public:
SumCountAccumulator(double initSum, size_t initCount): valueSum(initSum), sampleCount(initCount) {}
void addValue(double value)
{
if (value == HUGE_VAL || value != value)
{
return;
}
valueSum += value;
++sampleCount;
}
double valueSum;
size_t sampleCount;
};
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigFemNativeVisibleCellsStatCalc::valueSumAndSampleCount(size_t timeStepIndex, double& valueSum, size_t& sampleCount)
{
SumCountAccumulator acc(valueSum, sampleCount);
traverseElementNodes(acc, timeStepIndex);
valueSum = acc.valueSum;
sampleCount = acc.sampleCount;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigFemNativeVisibleCellsStatCalc::addDataToHistogramCalculator(size_t timeStepIndex, RigHistogramCalculator& histogramCalculator)
{
traverseElementNodes(histogramCalculator, timeStepIndex);
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
size_t RigFemNativeVisibleCellsStatCalc::timeStepCount()
{
return m_resultsData->frameCount();
}

View File

@@ -0,0 +1,106 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2015- Statoil ASA
// Copyright (C) 2015- Ceetron Solutions AS
//
// ResInsight 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.
//
// ResInsight 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 at <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#pragma once
//==================================================================================================
///
//==================================================================================================
#include "RigStatisticsCalculator.h"
#include "RigFemResultAddress.h"
#include "cvfArray.h"
class RigGeoMechCaseData;
class RigFemPartResultsCollection;
class RigFemNativeVisibleCellsStatCalc : public RigStatisticsCalculator
{
public:
RigFemNativeVisibleCellsStatCalc(RigGeoMechCaseData* femCase,
const RigFemResultAddress& resVarAddr,
const cvf::UByteArray* cellVisibilities);
virtual void minMaxCellScalarValues(size_t timeStepIndex, double& min, double& max);
virtual void posNegClosestToZero(size_t timeStepIndex, double& pos, double& neg);
virtual void valueSumAndSampleCount(size_t timeStepIndex, double& valueSum, size_t& sampleCount);
virtual void addDataToHistogramCalculator(size_t timeStepIndex, RigHistogramCalculator& histogramCalculator);
virtual size_t timeStepCount();
private:
RigGeoMechCaseData* m_caseData;
RigFemPartResultsCollection* m_resultsData;
RigFemResultAddress m_resVarAddr;
cvf::cref<cvf::UByteArray> m_cellVisibilities;
template <typename StatisticsAccumulator>
void traverseElementNodes(StatisticsAccumulator& accumulator, size_t timeStepIndex)
{
int partCount = m_caseData->femParts()->partCount();
if (m_resVarAddr.resultPosType == RIG_NODAL)
{
for (int pIdx = 0; pIdx < partCount; ++pIdx)
{
RigFemPart* part = m_caseData->femParts()->part(pIdx);
const std::vector<float>& values = m_resultsData->resultValues(m_resVarAddr, pIdx, (int)timeStepIndex);
int elmCount = part->elementCount();
for (int elmIdx = 0; elmIdx < elmCount; ++elmIdx)
{
if (!(*m_cellVisibilities)[elmIdx]) continue;
int elmNodeCount = RigFemTypes::elmentNodeCount(part->elementType(elmIdx));
for (int elmLocIdx = 0; elmLocIdx < elmNodeCount; ++elmLocIdx)
{
size_t elmNodeResIdx = part->elementNodeResultIdx(elmIdx, elmLocIdx);
int nodeIdx = part->nodeIdxFromElementNodeResultIdx(elmNodeResIdx);
accumulator.addValue(values[nodeIdx]);
}
}
}
}
else
{
for (int pIdx = 0; pIdx < partCount; ++pIdx)
{
RigFemPart* part = m_caseData->femParts()->part(pIdx);
const std::vector<float>& values = m_resultsData->resultValues(m_resVarAddr, pIdx, (int)timeStepIndex);
int elmCount = part->elementCount();
for (int elmIdx = 0; elmIdx < elmCount; ++elmIdx)
{
if (!(*m_cellVisibilities)[elmIdx]) continue;
int elmNodeCount = RigFemTypes::elmentNodeCount(part->elementType(elmIdx));
for (int elmLocIdx = 0; elmLocIdx < elmNodeCount; ++elmLocIdx)
{
size_t elmNodeResIdx = part->elementNodeResultIdx(elmIdx, elmLocIdx);
accumulator.addValue(values[elmNodeResIdx]);
}
}
}
}
}
};

View File

@@ -821,3 +821,11 @@ const std::vector<size_t>& RigFemPartResultsCollection::scalarValuesHistogram(co
return this->statistics(resVarAddr)->cellScalarValuesHistogram(frameIndex);
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
int RigFemPartResultsCollection::partCount() const
{
return m_femParts->partCount();
}

View File

@@ -43,7 +43,7 @@ public:
std::vector<std::string> stepNames();
bool assertResultsLoaded(const RigFemResultAddress& resVarAddr);
const std::vector<float>& resultValues(const RigFemResultAddress& resVarAddr, int partIndex, int frameIndex);
int partCount() const;
int frameCount();
@@ -67,7 +67,6 @@ private:
RigFemScalarResultFrames* calculateBarConvertedResult(int partIndex, const RigFemResultAddress &convertedResultAddr, const std::string fieldNameToConvert);
RigFemScalarResultFrames* calculateEnIpPorBarResult(int partIndex, const RigFemResultAddress &convertedResultAddr);
friend class RigFemNativeStatCalc;
cvf::Collection<RigFemPartResults> m_femPartResults;
cvf::ref<RifGeoMechReaderInterface> m_readerInterface;
cvf::cref<RigFemPartCollection> m_femParts;

View File

@@ -42,6 +42,7 @@
#include "RigFemPartResultsCollection.h"
#include "RigStatisticsDataCache.h"
#include "RigFemNativeVisibleCellsStatCalc.h"
CAF_PDM_SOURCE_INIT(Rim3dOverlayInfoConfig, "View3dOverlayInfoConfig");
//--------------------------------------------------------------------------------------------------
@@ -89,7 +90,7 @@ Rim3dOverlayInfoConfig::Rim3dOverlayInfoConfig()
CAF_PDM_InitFieldNoDefault(&m_statisticsTimeRange, "StatisticsTimeRange", "Statistics Time Range", "", "", "");
CAF_PDM_InitFieldNoDefault(&m_statisticsCellRange, "StatisticsCellRange", "Statistics Cell Range", "", "", "");
m_statisticsCellRange.uiCapability()->setUiHidden(true);
//m_statisticsCellRange.uiCapability()->setUiHidden(true);
}
//--------------------------------------------------------------------------------------------------
@@ -143,6 +144,8 @@ void Rim3dOverlayInfoConfig::update3DInfo()
m_viewDef->viewer()->showInfoText(showInfoText());
m_viewDef->viewer()->showHistogram(false);
m_viewDef->viewer()->showAnimationProgress(showAnimProgress());
m_isVisCellStatUpToDate = false;
RimEclipseView * reservoirView = dynamic_cast<RimEclipseView*>(m_viewDef.p());
if (reservoirView) updateEclipse3DInfo(reservoirView);
@@ -217,7 +220,8 @@ void Rim3dOverlayInfoConfig::updateEclipse3DInfo(RimEclipseView * reservoirView)
if (reservoirView->hasUserRequestedAnimation() && reservoirView->cellResult()->hasResult())
{
infoText += QString("<b>Cell Property:</b> %1 ").arg(propName);
// Wait until regression tests confirm new statisticks is ok infoText += QString("<br>Statistics for: ") + m_statisticsTimeRange().uiText() + " and " + m_statisticsCellRange().uiText();
// Wait until regression tests confirm new statisticks is ok :
//infoText += QString("<br>Statistics for: ") + m_statisticsTimeRange().uiText() + " and " + m_statisticsCellRange().uiText();
if (m_statisticsCellRange == ALL_CELLS)
{
@@ -349,116 +353,30 @@ void Rim3dOverlayInfoConfig::updateEclipse3DInfo(RimEclipseView * reservoirView)
//--------------------------------------------------------------------------------------------------
void Rim3dOverlayInfoConfig::updateGeoMech3DInfo(RimGeoMechView * geoMechView)
{
if (showInfoText())
RimGeoMechCase* geoMechCase = geoMechView->geoMechCase();
RigGeoMechCaseData* caseData = geoMechCase ? geoMechCase->geoMechData() : NULL;
bool isResultsInfoRelevant = caseData && geoMechView->hasUserRequestedAnimation() && geoMechView->cellResult()->hasResult();
// Retreive result stats if needed
double min = HUGE_VAL, max = HUGE_VAL;
double p10 = HUGE_VAL, p90 = HUGE_VAL;
double mean = HUGE_VAL;
const std::vector<size_t>* histogram = NULL;
if (showInfoText() || showHistogram())
{
QString infoText;
RimGeoMechCase* geoMechCase = geoMechView->geoMechCase();
RigGeoMechCaseData* caseData = geoMechCase ? geoMechCase->geoMechData() : NULL;
RigFemPartCollection* femParts = caseData ? caseData->femParts() : NULL;
if (femParts)
if (isResultsInfoRelevant)
{
QString caseName = geoMechCase->caseUserDescription();
QString cellCount = QString("%1").arg(femParts->totalElementCount());
QString zScale = QString::number(geoMechView->scaleZ());
infoText = QString(
"<p><b><center>-- %1 --</center></b><p>"
"<b>Cell count:</b> %2 <b>Z-Scale:</b> %3<br>").arg(caseName, cellCount, zScale);
if (geoMechView->hasUserRequestedAnimation() && geoMechView->cellResult()->hasResult())
RigFemResultAddress resAddress = geoMechView->cellResult()->resultAddress();
if (m_statisticsCellRange == ALL_CELLS)
{
QString resultPos;
QString fieldName = geoMechView->cellResult()->resultFieldUiName();
QString compName = geoMechView->cellResult()->resultComponentUiName();
if (!fieldName.isEmpty())
{
switch (geoMechView->cellResult()->resultPositionType())
{
case RIG_NODAL:
resultPos = "Nodal";
break;
case RIG_ELEMENT_NODAL:
resultPos = "Element nodal";
break;
case RIG_INTEGRATION_POINT:
resultPos = "Integration point";
break;
default:
break;
}
infoText += QString("<b>Cell result:</b> %1, %2, %3").arg(resultPos).arg(fieldName).arg(compName);
double min = HUGE_VAL, max = HUGE_VAL;
double p10 = HUGE_VAL, p90 = HUGE_VAL;
double mean = HUGE_VAL;
RigFemResultAddress resAddress = geoMechView->cellResult()->resultAddress();
if (m_statisticsTimeRange == ALL_TIMESTEPS)
{
caseData->femPartResults()->meanScalarValue (resAddress, &mean);
caseData->femPartResults()->minMaxScalarValues(resAddress,&min, &max);
caseData->femPartResults()->p10p90ScalarValues(resAddress, &p10, &p90);
}
else if (m_statisticsTimeRange == CURRENT_TIMESTEP)
{
int timeStepIdx = geoMechView->currentTimeStep();
caseData->femPartResults()->meanScalarValue (resAddress, timeStepIdx, &mean);
caseData->femPartResults()->minMaxScalarValues(resAddress, timeStepIdx, &min, &max);
caseData->femPartResults()->p10p90ScalarValues(resAddress, timeStepIdx, &p10, &p90);
}
infoText += QString("<table border=0 cellspacing=5 >"
"<tr> <td>Min</td> <td>P10</td> <td>Mean</td> <td>P90</td> <td>Max</td> </tr>"
"<tr> <td>%1</td> <td> %2</td> <td> %3</td> <td> %4</td> <td> %5</td> </tr>"
"</table>").arg(min).arg(p10).arg(mean).arg(p90).arg(max);
}
else
{
infoText += QString("<br>");
}
int currentTimeStep = geoMechView->currentTimeStep();
QString stepName = QString::fromStdString(caseData->femPartResults()->stepNames()[currentTimeStep]);
infoText += QString("<b>Time Step:</b> %1 <b>Time:</b> %2").arg(currentTimeStep).arg(stepName);
}
}
geoMechView->viewer()->setInfoText(infoText);
}
if (showHistogram())
{
if (geoMechView->hasUserRequestedAnimation() && geoMechView->cellResult()->hasResult())
{
geoMechView->viewer()->showHistogram(true);
// ToDo: Implement statistics for geomech data
RimGeoMechCase* geoMechCase = geoMechView->geoMechCase();
RigGeoMechCaseData* caseData = geoMechCase ? geoMechCase->geoMechData() : NULL;
if (caseData)
{
double min = HUGE_VAL, max = HUGE_VAL;
double p10 = HUGE_VAL, p90 = HUGE_VAL;
double mean = HUGE_VAL;
const std::vector<size_t>* histogram = NULL;
RigFemResultAddress resAddress = geoMechView->cellResult()->resultAddress();
if (m_statisticsTimeRange == ALL_TIMESTEPS)
{
caseData->femPartResults()->meanScalarValue(resAddress, &mean);
caseData->femPartResults()->minMaxScalarValues(resAddress, &min, &max);
caseData->femPartResults()->p10p90ScalarValues(resAddress, &p10, &p90);
histogram = &(caseData->femPartResults()->scalarValuesHistogram(resAddress));
}
else if (m_statisticsTimeRange == CURRENT_TIMESTEP)
@@ -467,16 +385,144 @@ void Rim3dOverlayInfoConfig::updateGeoMech3DInfo(RimGeoMechView * geoMechView)
caseData->femPartResults()->meanScalarValue(resAddress, timeStepIdx, &mean);
caseData->femPartResults()->minMaxScalarValues(resAddress, timeStepIdx, &min, &max);
caseData->femPartResults()->p10p90ScalarValues(resAddress, timeStepIdx, &p10, &p90);
histogram = &(caseData->femPartResults()->scalarValuesHistogram(resAddress, timeStepIdx));
}
else
{
CVF_ASSERT(false);
}
}
else if (m_statisticsCellRange == VISIBLE_CELLS)
{
this->updateVisCellStatsIfNeeded();
geoMechView->viewer()->setHistogram(min, max, *histogram);
geoMechView->viewer()->setHistogramPercentiles(p10, p90, mean);
if (m_statisticsTimeRange == ALL_TIMESTEPS)
{
// TODO: Only valid if we have no dynamic property filter
m_visibleCellStatistics->meanCellScalarValues(mean);
m_visibleCellStatistics->minMaxCellScalarValues(min, max);
m_visibleCellStatistics->p10p90CellScalarValues(p10, p90);
histogram = &(m_visibleCellStatistics->cellScalarValuesHistogram());
}
else if (m_statisticsTimeRange == CURRENT_TIMESTEP)
{
int timeStepIdx = geoMechView->currentTimeStep();
m_visibleCellStatistics->meanCellScalarValues(timeStepIdx, mean);
m_visibleCellStatistics->minMaxCellScalarValues(timeStepIdx, min, max);
m_visibleCellStatistics->p10p90CellScalarValues(timeStepIdx, p10, p90);
histogram = &(m_visibleCellStatistics->cellScalarValuesHistogram(timeStepIdx));
}
}
}
}
// Compose text
if (showInfoText())
{
QString infoText;
RigFemPartCollection* femParts = caseData ? caseData->femParts() : NULL;
if (femParts)
{
QString caseName = geoMechCase->caseUserDescription();
QString cellCount = QString("%1").arg(femParts->totalElementCount());
QString zScale = QString::number(geoMechView->scaleZ());
infoText = QString(
"<p><b><center>-- %1 --</center></b><p>"
"<b>Cell count:</b> %2 <b>Z-Scale:</b> %3<br>").arg(caseName, cellCount, zScale);
}
if (isResultsInfoRelevant)
{
{
QString resultPos;
QString fieldName = geoMechView->cellResult()->resultFieldUiName();
QString compName = geoMechView->cellResult()->resultComponentUiName();
switch (geoMechView->cellResult()->resultPositionType())
{
case RIG_NODAL:
resultPos = "Nodal";
break;
case RIG_ELEMENT_NODAL:
resultPos = "Element nodal";
break;
case RIG_INTEGRATION_POINT:
resultPos = "Integration point";
break;
default:
break;
}
infoText += QString("<b>Cell result:</b> %1, %2, %3").arg(resultPos).arg(fieldName).arg(compName);
}
{
infoText += QString("<table border=0 cellspacing=5 >"
"<tr> <td>Min</td> <td>P10</td> <td>Mean</td> <td>P90</td> <td>Max</td> </tr>"
"<tr> <td>%1</td> <td> %2</td> <td> %3</td> <td> %4</td> <td> %5</td> </tr>"
"</table>").arg(min).arg(p10).arg(mean).arg(p90).arg(max);
}
{
int currentTimeStep = geoMechView->currentTimeStep();
QString stepName = QString::fromStdString(caseData->femPartResults()->stepNames()[currentTimeStep]);
infoText += QString("<b>Time Step:</b> %1 <b>Time:</b> %2").arg(currentTimeStep).arg(stepName);
}
}
geoMechView->viewer()->setInfoText(infoText);
}
// Populate histogram
if (showHistogram())
{
if (isResultsInfoRelevant)
{
geoMechView->viewer()->showHistogram(true);
geoMechView->viewer()->setHistogram(min, max, *histogram);
geoMechView->viewer()->setHistogramPercentiles(p10, p90, mean);
}
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void Rim3dOverlayInfoConfig::updateVisCellStatsIfNeeded()
{
RimEclipseView * eclipseView = dynamic_cast<RimEclipseView*>(m_viewDef.p());
RimGeoMechView * geoMechView = dynamic_cast<RimGeoMechView*>(m_viewDef.p());
if (!m_isVisCellStatUpToDate)
{
cvf::ref<RigStatisticsCalculator> calc;
if (geoMechView)
{
RigFemResultAddress resAddress = geoMechView->cellResult()->resultAddress();
calc = new RigFemNativeVisibleCellsStatCalc(geoMechView->geoMechCase()->geoMechData(),
resAddress,
geoMechView->currentTotalCellVisibility().p());
m_visibleCellStatistics = new RigStatisticsDataCache(calc.p());
m_isVisCellStatUpToDate = true;
}
else if (eclipseView)
{
// RigFemResultAddress resAddress = geoMechView->cellResult()->resultAddress();
// cvf::ref<RigEclipseNativeVisibleCellsStatCalc> calc = new RigEclipseNativeVisibleCellsStatCalc(geoMechView->geoMechCase()->geoMechData(),
// resAddress,
// geoMechView->currentTotalCellVisibility().p());
}
m_visibleCellStatistics = new RigStatisticsDataCache(calc.p());
m_isVisCellStatUpToDate = true;
}
}

View File

@@ -81,6 +81,9 @@ private:
cvf::Vec2ui m_position;
void updateVisCellStatsIfNeeded();
bool m_isVisCellStatUpToDate;
cvf::ref<RigStatisticsDataCache> m_visibleCellStatistics;
};

View File

@@ -31,6 +31,8 @@
RigStatisticsDataCache::RigStatisticsDataCache(RigStatisticsCalculator* statisticsCalculator)
: m_statisticsCalculator(statisticsCalculator)
{
CVF_ASSERT(m_statisticsCalculator.notNull());
clearAllStatistics();
}
@@ -153,6 +155,40 @@ void RigStatisticsDataCache::posNegClosestToZero(size_t timeStepIndex, double& p
posNearZero = m_statsPrTs[timeStepIndex].m_posClosestToZero;
negNearZero = m_statsPrTs[timeStepIndex].m_negClosestToZero;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigStatisticsDataCache::meanCellScalarValues(double& meanValue)
{
if (!m_statsAllTimesteps.m_isMeanCalculated)
{
m_statisticsCalculator->meanCellScalarValue(m_statsAllTimesteps.m_meanValue);
m_statsAllTimesteps.m_isMeanCalculated = true;
}
meanValue = m_statsAllTimesteps.m_meanValue;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigStatisticsDataCache::meanCellScalarValues(size_t timeStepIndex, double& meanValue)
{
if (timeStepIndex >= m_statsPrTs.size())
{
m_statsPrTs.resize(timeStepIndex + 1);
}
if (!m_statsPrTs[timeStepIndex].m_isMeanCalculated)
{
m_statisticsCalculator->meanCellScalarValue(timeStepIndex, m_statsPrTs[timeStepIndex].m_meanValue);
m_statsPrTs[timeStepIndex].m_isMeanCalculated = true;
}
meanValue = m_statsPrTs[timeStepIndex].m_meanValue;
}
//--------------------------------------------------------------------------------------------------
///
@@ -196,34 +232,6 @@ void RigStatisticsDataCache::p10p90CellScalarValues(size_t timeStepIndex, double
p90 = m_statsPrTs[timeStepIndex].m_p90;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigStatisticsDataCache::meanCellScalarValues(double& meanValue)
{
if (!m_statsAllTimesteps.m_isMeanCalculated)
{
m_statisticsCalculator->meanCellScalarValue(m_statsAllTimesteps.m_meanValue);
m_statsAllTimesteps.m_isMeanCalculated = true;
}
meanValue = m_statsAllTimesteps.m_meanValue;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigStatisticsDataCache::meanCellScalarValues(size_t timeStepIndex, double& meanValue)
{
if (!m_statsPrTs[timeStepIndex].m_isMeanCalculated)
{
m_statisticsCalculator->meanCellScalarValue(timeStepIndex, m_statsPrTs[timeStepIndex].m_meanValue);
m_statsPrTs[timeStepIndex].m_isMeanCalculated = true;
}
meanValue = m_statsPrTs[timeStepIndex].m_meanValue;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------

View File

@@ -186,7 +186,28 @@ RigHistogramCalculator::RigHistogramCalculator(double min, double max, size_t nB
for (size_t i = 0; i < m_histogram->size(); ++i) (*m_histogram)[i] = 0;
m_range = max - min;
maxIndex = nBins-1;
m_maxIndex = nBins-1;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigHistogramCalculator::addValue(double value)
{
if (value == HUGE_VAL || value != value)
{
return;
}
size_t index = 0;
if (m_maxIndex > 0) index = (size_t)(m_maxIndex*(value - m_min)/m_range);
if (index < m_histogram->size()) // Just clip to the max min range (-index will overflow to positive )
{
(*m_histogram)[index]++;
m_observationCount++;
}
}
//--------------------------------------------------------------------------------------------------
@@ -197,20 +218,7 @@ void RigHistogramCalculator::addData(const std::vector<double>& data)
assert(m_histogram);
for (size_t i = 0; i < data.size(); ++i)
{
if (data[i] == HUGE_VAL)
{
continue;
}
size_t index = 0;
if (maxIndex > 0) index = (size_t)(maxIndex*(data[i] - m_min)/m_range);
if(index < m_histogram->size()) // Just clip to the max min range (-index will overflow to positive )
{
(*m_histogram)[index]++;
m_observationCount++;
}
addValue(data[i]);
}
}
@@ -222,20 +230,7 @@ void RigHistogramCalculator::addData(const std::vector<float>& data)
assert(m_histogram);
for (size_t i = 0; i < data.size(); ++i)
{
if (data[i] == HUGE_VAL)
{
continue;
}
size_t index = 0;
if (maxIndex > 0) index = (size_t)(maxIndex*(data[i] - m_min)/m_range);
if(index < m_histogram->size()) // Just clip to the max min range (-index will overflow to positive )
{
(*m_histogram)[index]++;
m_observationCount++;
}
addValue(data[i]);
}
}

View File

@@ -41,6 +41,8 @@ public:
void addData(const std::vector<double>& data);
void addData(const std::vector<float>& data);
void addValue(double value);
/// Calculates the estimated percentile from the histogram.
/// the percentile is the domain value at which pVal of the observations are below it.
/// Will only consider observed values between min and max, as all other values are discarded from the histogram
@@ -48,7 +50,7 @@ public:
double calculatePercentil(double pVal);
private:
size_t maxIndex;
size_t m_maxIndex;
double m_range;
double m_min;
size_t m_observationCount;