OilfieldToolsNet/ResInsight
4
0
Fork 0
mirror of https://github.com/OPM/ResInsight.git synced 2025-02-25 18:55:39 -06:00
Code Issues Packages Projects Releases Wiki Activity

#2659 Ensemble statistics. Move classes to separate files

Browse Source
This commit is contained in:
Bjørn Erik Jensen 2018-06-07 10:50:32 +02:00
parent 12aec1eba2
commit cb6657266a
8 changed files with 414 additions and 308 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
ApplicationCode/FileInterface/CMakeLists_files.cmake
View File

@ -42,6 +42,7 @@ ${CMAKE_CURRENT_LIST_DIR}/RifSummaryCaseRestartSelector.h
${CMAKE_CURRENT_LIST_DIR}/RifEnsembleParametersReader.h
${CMAKE_CURRENT_LIST_DIR}/RifCaseRealizationParametersReader.h
${CMAKE_CURRENT_LIST_DIR}/RifFileParseTools.h
${CMAKE_CURRENT_LIST_DIR}/RifEnsembleStatisticsReader.h
# HDF5 file reader is directly included in ResInsight main CmakeList.txt
#${CMAKE_CURRENT_LIST_DIR}/RifHdf5Reader.h
@ -88,6 +89,7 @@ ${CMAKE_CURRENT_LIST_DIR}/RifSummaryCaseRestartSelector.cpp
${CMAKE_CURRENT_LIST_DIR}/RifEnsembleParametersReader.cpp
${CMAKE_CURRENT_LIST_DIR}/RifCaseRealizationParametersReader.cpp
${CMAKE_CURRENT_LIST_DIR}/RifFileParseTools.cpp
${CMAKE_CURRENT_LIST_DIR}/RifEnsembleStatisticsReader.cpp
# HDF5 file reader is directly included in ResInsight main CmakeList.txt
#${CMAKE_CURRENT_LIST_DIR}/RifHdf5Reader.cpp

86
ApplicationCode/FileInterface/RifEnsembleStatisticsReader.cpp Normal file
View File

@ -0,0 +1,86 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2017- Statoil ASA
//
// 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 "RifEnsembleStatisticsReader.h"
#include "RimEnsembleStatisticsCase.h"
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RifEnsembleStatisticsReader::RifEnsembleStatisticsReader(RimEnsembleStatisticsCase* ensStatCase)
{
CVF_ASSERT(ensStatCase);
m_ensembleStatCase = ensStatCase;
m_allResultAddresses = std::vector<RifEclipseSummaryAddress>(
{
RifEclipseSummaryAddress::ensembleStatisticsAddress(ENSEMBLE_STAT_P10_QUANTITY_NAME),
RifEclipseSummaryAddress::ensembleStatisticsAddress(ENSEMBLE_STAT_P50_QUANTITY_NAME),
RifEclipseSummaryAddress::ensembleStatisticsAddress(ENSEMBLE_STAT_P90_QUANTITY_NAME),
RifEclipseSummaryAddress::ensembleStatisticsAddress(ENSEMBLE_STAT_MEAN_QUANTITY_NAME)
});
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<time_t>& RifEnsembleStatisticsReader::timeSteps(const RifEclipseSummaryAddress& resultAddress) const
{
return m_ensembleStatCase->timeSteps();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
bool RifEnsembleStatisticsReader::values(const RifEclipseSummaryAddress& resultAddress, std::vector<double>* values) const
{
if (!validateAddress(resultAddress)) return false;
const std::vector<double>* sourceData = nullptr;
if (resultAddress.quantityName() == ENSEMBLE_STAT_P10_QUANTITY_NAME) sourceData = &m_ensembleStatCase->p10();
else if (resultAddress.quantityName() == ENSEMBLE_STAT_P50_QUANTITY_NAME) sourceData = &m_ensembleStatCase->p50();
else if (resultAddress.quantityName() == ENSEMBLE_STAT_P90_QUANTITY_NAME) sourceData = &m_ensembleStatCase->p90();
else if (resultAddress.quantityName() == ENSEMBLE_STAT_MEAN_QUANTITY_NAME) sourceData = &m_ensembleStatCase->mean();
if (!sourceData) return false;
values->clear();
values->reserve(sourceData->size());
for (auto val : *sourceData) values->push_back(val);
return true;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::string RifEnsembleStatisticsReader::unitName(const RifEclipseSummaryAddress& resultAddress) const
{
return "(RifEnsembleStatisticsReader::unitName)";
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
bool RifEnsembleStatisticsReader::validateAddress(const RifEclipseSummaryAddress& address) const
{
return address.category() == RifEclipseSummaryAddress::SUMMARY_ENSEMBLE_STATISTICS &&
!address.quantityName().empty();
}

45
ApplicationCode/FileInterface/RifEnsembleStatisticsReader.h Normal file
View File

@ -0,0 +1,45 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2017- Statoil ASA
//
// 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 "RifEclipseSummaryAddress.h"
#include "RifSummaryReaderInterface.h"
class RimEnsembleStatisticsCase;
//==================================================================================================
///
//==================================================================================================
class RifEnsembleStatisticsReader : public RifSummaryReaderInterface
{
public:
RifEnsembleStatisticsReader(RimEnsembleStatisticsCase* ensStatCase);
virtual const std::vector<time_t>& timeSteps(const RifEclipseSummaryAddress& resultAddress) const override;
virtual bool values(const RifEclipseSummaryAddress& resultAddress, std::vector<double>* values) const override;
virtual std::string unitName(const RifEclipseSummaryAddress& resultAddress) const override;
private:
bool validateAddress(const RifEclipseSummaryAddress& address) const;
private:
RimEnsembleStatisticsCase * m_ensembleStatCase;
};

2
ApplicationCode/ProjectDataModel/Summary/CMakeLists_files.cmake
View File

@ -35,6 +35,7 @@ ${CMAKE_CURRENT_LIST_DIR}/RimEnsembleCurveSetColorManager.h
${CMAKE_CURRENT_LIST_DIR}/RimEnsembleCurveFilter.h
${CMAKE_CURRENT_LIST_DIR}/RimEnsembleCurveFilterCollection.h
${CMAKE_CURRENT_LIST_DIR}/RimEnsembleStatistics.h
${CMAKE_CURRENT_LIST_DIR}/RimEnsembleStatisticsCase.h
)
set (SOURCE_GROUP_SOURCE_FILES
@ -73,6 +74,7 @@ ${CMAKE_CURRENT_LIST_DIR}/RimEnsembleCurveSetColorManager.cpp
${CMAKE_CURRENT_LIST_DIR}/RimEnsembleCurveFilter.cpp
${CMAKE_CURRENT_LIST_DIR}/RimEnsembleCurveFilterCollection.cpp
${CMAKE_CURRENT_LIST_DIR}/RimEnsembleStatistics.cpp
${CMAKE_CURRENT_LIST_DIR}/RimEnsembleStatisticsCase.cpp
)
list(APPEND CODE_HEADER_FILES

255
ApplicationCode/ProjectDataModel/Summary/RimEnsembleCurveSet.cpp
View File

@ -23,6 +23,7 @@
#include "SummaryPlotCommands/RicSummaryCurveCreator.h"
#include "RifEnsembleStatisticsReader.h"
#include "RifReaderEclipseSummary.h"
#include "RigStatisticsMath.h"
@ -33,6 +34,7 @@
#include "RimEnsembleCurveSetCollection.h"
#include "RimEnsembleCurveSetColorManager.h"
#include "RimEnsembleStatistics.h"
#include "RimEnsembleStatisticsCase.h"
#include "RimProject.h"
#include "RimRegularLegendConfig.h"
#include "RimSummaryAddress.h"
@ -302,8 +304,8 @@ std::vector<RimSummaryCurve*> RimEnsembleCurveSet::curves() const
//--------------------------------------------------------------------------------------------------
void RimEnsembleCurveSet::deleteEnsembleCurves()
{
std::vector<int> curvesIndexesToDelete;
for (int c = 0; c < m_curves.size(); c++)
std::vector<size_t> curvesIndexesToDelete;
for (size_t c = 0; c < m_curves.size(); c++)
{
RimSummaryCurve* curve = m_curves[c];
if (curve->summaryAddressY().category() != RifEclipseSummaryAddress::SUMMARY_ENSEMBLE_STATISTICS)
@ -312,7 +314,7 @@ void RimEnsembleCurveSet::deleteEnsembleCurves()
while (curvesIndexesToDelete.size() > 0)
{
int currIndex = curvesIndexesToDelete.back();
size_t currIndex = curvesIndexesToDelete.back();
delete m_curves[currIndex];
m_curves.erase(currIndex);
curvesIndexesToDelete.pop_back();
@ -324,8 +326,8 @@ void RimEnsembleCurveSet::deleteEnsembleCurves()
//--------------------------------------------------------------------------------------------------
void RimEnsembleCurveSet::deleteStatisticsCurves()
{
std::vector<int> curvesIndexesToDelete;
for (int c = 0; c < m_curves.size(); c++)
std::vector<size_t> curvesIndexesToDelete;
for (size_t c = 0; c < m_curves.size(); c++)
{
RimSummaryCurve* curve = m_curves[c];
if (curve->summaryAddressY().category() == RifEclipseSummaryAddress::SUMMARY_ENSEMBLE_STATISTICS)
@ -334,7 +336,7 @@ void RimEnsembleCurveSet::deleteStatisticsCurves()
while (curvesIndexesToDelete.size() > 0)
{
int currIndex = curvesIndexesToDelete.back();
size_t currIndex = curvesIndexesToDelete.back();
delete m_curves[currIndex];
m_curves.erase(currIndex);
curvesIndexesToDelete.pop_back();
@ -1139,244 +1141,3 @@ void RimEnsembleCurveSet::updateLegendMappingMode()
break;
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RifEnsembleStatisticsReader::RifEnsembleStatisticsReader(RimEnsembleStatisticsCase* ensStatCase)
{
CVF_ASSERT(ensStatCase);
m_ensembleStatCase = ensStatCase;
m_allResultAddresses = std::vector<RifEclipseSummaryAddress>(
{
RifEclipseSummaryAddress::ensembleStatisticsAddress(ENSEMBLE_STAT_P10_QUANTITY_NAME),
RifEclipseSummaryAddress::ensembleStatisticsAddress(ENSEMBLE_STAT_P50_QUANTITY_NAME),
RifEclipseSummaryAddress::ensembleStatisticsAddress(ENSEMBLE_STAT_P90_QUANTITY_NAME),
RifEclipseSummaryAddress::ensembleStatisticsAddress(ENSEMBLE_STAT_MEAN_QUANTITY_NAME)
});
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<time_t>& RifEnsembleStatisticsReader::timeSteps(const RifEclipseSummaryAddress& resultAddress) const
{
return m_ensembleStatCase->timeSteps();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
bool RifEnsembleStatisticsReader::values(const RifEclipseSummaryAddress& resultAddress, std::vector<double>* values) const
{
if (!validateAddress(resultAddress)) return false;
const std::vector<double>* sourceData = nullptr;
if (resultAddress.quantityName() == ENSEMBLE_STAT_P10_QUANTITY_NAME) sourceData = &m_ensembleStatCase->p10();
else if (resultAddress.quantityName() == ENSEMBLE_STAT_P50_QUANTITY_NAME) sourceData = &m_ensembleStatCase->p50();
else if (resultAddress.quantityName() == ENSEMBLE_STAT_P90_QUANTITY_NAME) sourceData = &m_ensembleStatCase->p90();
else if (resultAddress.quantityName() == ENSEMBLE_STAT_MEAN_QUANTITY_NAME) sourceData = &m_ensembleStatCase->mean();
if (!sourceData) return false;
values->clear();
values->reserve(sourceData->size());
for (auto val : *sourceData) values->push_back(val);
return true;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::string RifEnsembleStatisticsReader::unitName(const RifEclipseSummaryAddress& resultAddress) const
{
return "(RifEnsembleStatisticsReader::unitName)";
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
bool RifEnsembleStatisticsReader::validateAddress(const RifEclipseSummaryAddress& address) const
{
return address.category() == RifEclipseSummaryAddress::SUMMARY_ENSEMBLE_STATISTICS &&
!address.quantityName().empty();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RimEnsembleStatisticsCase::RimEnsembleStatisticsCase(RimEnsembleCurveSet* curveSet)
{
m_curveSet = curveSet;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<time_t>& RimEnsembleStatisticsCase::timeSteps() const
{
return m_timeSteps;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<double>& RimEnsembleStatisticsCase::p10() const
{
return m_p10Data;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<double>& RimEnsembleStatisticsCase::p50() const
{
return m_p50Data;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<double>& RimEnsembleStatisticsCase::p90() const
{
return m_p90Data;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<double>& RimEnsembleStatisticsCase::mean() const
{
return m_meanData;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
QString RimEnsembleStatisticsCase::caseName()
{
return "Ensemble Statistics";
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RimEnsembleStatisticsCase::createSummaryReaderInterface()
{
m_statisticsReader.reset(new RifEnsembleStatisticsReader(this));
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RifSummaryReaderInterface* RimEnsembleStatisticsCase::summaryReader()
{
return m_statisticsReader.get();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const RimEnsembleCurveSet* RimEnsembleStatisticsCase::curveSet() const
{
return m_curveSet;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RimEnsembleStatisticsCase::calculate(const RimSummaryCaseCollection* ensemble, const RifEclipseSummaryAddress& inputAddress)
{
RigTimeHistoryCurveMerger curveMerger;
int caseCount = (int)ensemble->allSummaryCases().size();
{
for (const auto& sumCase : ensemble->allSummaryCases())
{
const auto& reader = sumCase->summaryReader();
if (reader)
{
std::vector<time_t> timeSteps = reader->timeSteps(inputAddress);
std::vector<double> values;
reader->values(inputAddress, &values);
curveMerger.addCurveData(values, timeSteps);
}
}
curveMerger.computeInterpolatedValues();
}
const std::vector<time_t>& allTimeSteps = curveMerger.allTimeSteps();
std::vector<std::vector<double>> allValues;
{
for (int c = 0; c < caseCount; c++)
{
allValues.push_back(curveMerger.interpolatedCurveValuesForAllTimeSteps(c));
}
}
clearData();
m_timeSteps = allTimeSteps;
for (int t = 0; t < (int)allTimeSteps.size(); t++)
{
std::vector<double> valuesForTimeSteps;
valuesForTimeSteps.reserve(caseCount);
for (int c = 0; c < caseCount; c++)
{
valuesForTimeSteps.push_back(allValues[c][t]);
}
{
double min, max, range, mean, stdev;
RigStatisticsMath::calculateBasicStatistics(valuesForTimeSteps, &min, &max, nullptr, &range, &mean, &stdev);
std::vector<size_t> histogram;
RigHistogramCalculator histCalc(min, max, 100, &histogram);
histCalc.addData(valuesForTimeSteps);
double p10, p50, p90;
p10 = histCalc.calculatePercentil(0.1);
p50 = histCalc.calculatePercentil(0.5);
p90 = histCalc.calculatePercentil(0.9);
m_p10Data.push_back(p10);
m_p50Data.push_back(p50);
m_p90Data.push_back(p90);
m_meanData.push_back(mean);
}
}
m_addressUsedInLastCalculation = inputAddress;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RimEnsembleStatisticsCase::calculate()
{
auto inputAddress = m_curveSet->summaryAddress();
auto ensemble = m_curveSet->summaryCaseCollection();
if (m_statisticsReader && ensemble && inputAddress.isValid())
{
calculate(ensemble, inputAddress);
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RimEnsembleStatisticsCase::clearData()
{
m_timeSteps.clear();
m_p10Data.clear();
m_p50Data.clear();
m_p90Data.clear();
m_meanData.clear();
m_addressUsedInLastCalculation = RifEclipseSummaryAddress();
}

61
ApplicationCode/ProjectDataModel/Summary/RimEnsembleCurveSet.h
View File

@ -55,67 +55,6 @@ class QKeyEvent;
class RimEnsembleStatisticsCase;
//==================================================================================================
///
//==================================================================================================
class RifEnsembleStatisticsReader : public RifSummaryReaderInterface
{
public:
RifEnsembleStatisticsReader(RimEnsembleStatisticsCase* ensStatCase);
virtual const std::vector<time_t>& timeSteps(const RifEclipseSummaryAddress& resultAddress) const override;
virtual bool values(const RifEclipseSummaryAddress& resultAddress, std::vector<double>* values) const override;
virtual std::string unitName(const RifEclipseSummaryAddress& resultAddress) const override;
private:
bool validateAddress(const RifEclipseSummaryAddress& address) const;
private:
RimEnsembleStatisticsCase * m_ensembleStatCase;
};
//==================================================================================================
///
//==================================================================================================
class RimEnsembleStatisticsCase : public RimSummaryCase
{
public:
RimEnsembleStatisticsCase(RimEnsembleCurveSet* curveSet);
const std::vector<time_t>& timeSteps() const;
const std::vector<double>& p10() const;
const std::vector<double>& p50() const;
const std::vector<double>& p90() const;
const std::vector<double>& mean() const;
virtual QString caseName() override;
virtual void createSummaryReaderInterface() override;
virtual RifSummaryReaderInterface* summaryReader() override;
virtual void updateFilePathsFromProjectPath(const QString& newProjectPath, const QString& oldProjectPath) override {}
const RimEnsembleCurveSet* curveSet() const;
void calculate();
void calculate(const RimSummaryCaseCollection* ensemble, const RifEclipseSummaryAddress& inputAddress);
private:
void clearData();
private:
std::unique_ptr<RifEnsembleStatisticsReader> m_statisticsReader;
RimEnsembleCurveSet* m_curveSet;
std::vector<time_t> m_timeSteps;
std::vector<double> m_p10Data;
std::vector<double> m_p50Data;
std::vector<double> m_p90Data;
std::vector<double> m_meanData;
RifEclipseSummaryAddress m_addressUsedInLastCalculation;
};
//==================================================================================================
///
//==================================================================================================

202
ApplicationCode/ProjectDataModel/Summary/RimEnsembleStatisticsCase.cpp Normal file
View File

@ -0,0 +1,202 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2017- Statoil ASA
//
// 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 "RimEnsembleStatisticsCase.h"
#include "RifEnsembleStatisticsReader.h"
#include "RigStatisticsMath.h"
#include "RigTimeHistoryCurveMerger.h"
#include "RimEnsembleCurveSet.h"
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RimEnsembleStatisticsCase::RimEnsembleStatisticsCase(RimEnsembleCurveSet* curveSet)
{
m_curveSet = curveSet;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<time_t>& RimEnsembleStatisticsCase::timeSteps() const
{
return m_timeSteps;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<double>& RimEnsembleStatisticsCase::p10() const
{
return m_p10Data;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<double>& RimEnsembleStatisticsCase::p50() const
{
return m_p50Data;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<double>& RimEnsembleStatisticsCase::p90() const
{
return m_p90Data;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<double>& RimEnsembleStatisticsCase::mean() const
{
return m_meanData;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
QString RimEnsembleStatisticsCase::caseName()
{
return "Ensemble Statistics";
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RimEnsembleStatisticsCase::createSummaryReaderInterface()
{
m_statisticsReader.reset(new RifEnsembleStatisticsReader(this));
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RifSummaryReaderInterface* RimEnsembleStatisticsCase::summaryReader()
{
return m_statisticsReader.get();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const RimEnsembleCurveSet* RimEnsembleStatisticsCase::curveSet() const
{
return m_curveSet;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RimEnsembleStatisticsCase::calculate(const RimSummaryCaseCollection* ensemble, const RifEclipseSummaryAddress& inputAddress)
{
RigTimeHistoryCurveMerger curveMerger;
int caseCount = (int)ensemble->allSummaryCases().size();
{
for (const auto& sumCase : ensemble->allSummaryCases())
{
const auto& reader = sumCase->summaryReader();
if (reader)
{
std::vector<time_t> timeSteps = reader->timeSteps(inputAddress);
std::vector<double> values;
reader->values(inputAddress, &values);
curveMerger.addCurveData(values, timeSteps);
}
}
curveMerger.computeInterpolatedValues();
}
const std::vector<time_t>& allTimeSteps = curveMerger.allTimeSteps();
std::vector<std::vector<double>> allValues;
{
for (int c = 0; c < caseCount; c++)
{
allValues.push_back(curveMerger.interpolatedCurveValuesForAllTimeSteps(c));
}
}
clearData();
m_timeSteps = allTimeSteps;
for (int t = 0; t < (int)allTimeSteps.size(); t++)
{
std::vector<double> valuesForTimeSteps;
valuesForTimeSteps.reserve(caseCount);
for (int c = 0; c < caseCount; c++)
{
valuesForTimeSteps.push_back(allValues[c][t]);
}
{
double min, max, range, mean, stdev;
RigStatisticsMath::calculateBasicStatistics(valuesForTimeSteps, &min, &max, nullptr, &range, &mean, &stdev);
std::vector<size_t> histogram;
RigHistogramCalculator histCalc(min, max, 100, &histogram);
histCalc.addData(valuesForTimeSteps);
double p10, p50, p90;
p10 = histCalc.calculatePercentil(0.1);
p50 = histCalc.calculatePercentil(0.5);
p90 = histCalc.calculatePercentil(0.9);
m_p10Data.push_back(p10);
m_p50Data.push_back(p50);
m_p90Data.push_back(p90);
m_meanData.push_back(mean);
}
}
m_addressUsedInLastCalculation = inputAddress;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RimEnsembleStatisticsCase::calculate()
{
auto inputAddress = m_curveSet->summaryAddress();
auto ensemble = m_curveSet->summaryCaseCollection();
if (m_statisticsReader && ensemble && inputAddress.isValid())
{
calculate(ensemble, inputAddress);
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RimEnsembleStatisticsCase::clearData()
{
m_timeSteps.clear();
m_p10Data.clear();
m_p50Data.clear();
m_p90Data.clear();
m_meanData.clear();
m_addressUsedInLastCalculation = RifEclipseSummaryAddress();
}

69
ApplicationCode/ProjectDataModel/Summary/RimEnsembleStatisticsCase.h Normal file
View File

@ -0,0 +1,69 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2017- Statoil ASA
//
// 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 "RifEclipseSummaryAddress.h"
#include "RimSummaryCase.h"
class RifEnsembleStatisticsReader;
class RimEnsembleCurveSet;
//==================================================================================================
///
//==================================================================================================
class RimEnsembleStatisticsCase : public RimSummaryCase
{
public:
RimEnsembleStatisticsCase(RimEnsembleCurveSet* curveSet);
const std::vector<time_t>& timeSteps() const;
const std::vector<double>& p10() const;
const std::vector<double>& p50() const;
const std::vector<double>& p90() const;
const std::vector<double>& mean() const;
virtual QString caseName() override;
virtual void createSummaryReaderInterface() override;
virtual RifSummaryReaderInterface* summaryReader() override;
virtual void updateFilePathsFromProjectPath(const QString& newProjectPath, const QString& oldProjectPath) override {}
const RimEnsembleCurveSet* curveSet() const;
void calculate();
void calculate(const RimSummaryCaseCollection* ensemble, const RifEclipseSummaryAddress& inputAddress);
private:
void clearData();
private:
std::unique_ptr<RifEnsembleStatisticsReader> m_statisticsReader;
RimEnsembleCurveSet* m_curveSet;
std::vector<time_t> m_timeSteps;
std::vector<double> m_p10Data;
std::vector<double> m_p50Data;
std::vector<double> m_p90Data;
std::vector<double> m_meanData;
RifEclipseSummaryAddress m_addressUsedInLastCalculation;
};