#7493 Compute statistics for fracture group.

This commit is contained in:
Kristian Bendiksen 2021-03-17 19:04:04 +01:00
parent 6f9f77e009
commit e3852e5d6e
2 changed files with 370 additions and 0 deletions

View File

@ -18,8 +18,23 @@
#include "RimFractureGroupStatistics.h"
#include "RiaLogging.h"
#include "RiaPreferences.h"
#include "RigFractureGrid.h"
#include "RigSlice2D.h"
#include "RigStatisticsMath.h"
#include "RigStimPlanFractureDefinition.h"
#include "RifCsvDataTableFormatter.h"
#include "RifStimPlanXmlReader.h"
#include "cafPdmUiTextEditor.h"
#include <cmath>
#include <QFile>
CAF_PDM_SOURCE_INIT( RimFractureGroupStatistics, "FractureGroupStatistics" );
//--------------------------------------------------------------------------------------------------
@ -92,4 +107,320 @@ void RimFractureGroupStatistics::defineEditorAttribute( const caf::PdmFieldHandl
void RimFractureGroupStatistics::loadAndUpdateData()
{
m_filePathsTable = generateFilePathsTable();
computeStatistics();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RimFractureGroupStatistics::computeStatistics()
{
auto unitSystem = RiaDefines::EclipseUnitSystem::UNITS_METRIC;
std::vector<cvf::ref<RigStimPlanFractureDefinition>> stimPlanFractureDefinitions =
readFractureDefinitions( m_filePaths.v(), unitSystem );
std::set<QString> availableResults = findAllResultNames( stimPlanFractureDefinitions );
std::vector<cvf::cref<RigFractureGrid>> fractureGrids = createFractureGrids( stimPlanFractureDefinitions, unitSystem );
auto [minX, maxX, minY, maxY] = findExtentsOfGrids( fractureGrids );
int numSamplesX = 100;
int numSamplesY = 150;
double sampleDistanceX = ( maxX - minX ) / numSamplesX;
double sampleDistanceY = ( maxY - minY ) / numSamplesY;
RiaLogging::info(
QString( "Fracture Group Size: X = [%1, %2] Y = [%3, %4]" ).arg( minX ).arg( maxX ).arg( minY ).arg( maxY ) );
RiaLogging::info( QString( "Output size: %1x%2. Sampling Distance X = %3 Sampling Distance Y = %4" )
.arg( numSamplesX )
.arg( numSamplesY )
.arg( sampleDistanceX )
.arg( sampleDistanceY ) );
std::vector<std::vector<double>> samples( numSamplesX * numSamplesY );
sampleAllGrids( fractureGrids, samples, minX, minY, numSamplesX, numSamplesY, sampleDistanceX, sampleDistanceY );
std::map<QString, std::shared_ptr<RigSlice2D>> statisticsGrids;
generateStatisticsGrids( samples, numSamplesX, numSamplesY, statisticsGrids );
writeStatisticsToCsv( "/tmp/mean.csv", *statisticsGrids["mean"] );
writeStatisticsToCsv( "/tmp/min.csv", *statisticsGrids["min"] );
writeStatisticsToCsv( "/tmp/max.csv", *statisticsGrids["max"] );
writeStatisticsToCsv( "/tmp/p10.csv", *statisticsGrids["p10"] );
writeStatisticsToCsv( "/tmp/p50.csv", *statisticsGrids["p50"] );
writeStatisticsToCsv( "/tmp/p90.csv", *statisticsGrids["p90"] );
writeStatisticsToCsv( "/tmp/occurrence.csv", *statisticsGrids["occurrence"] );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::vector<cvf::ref<RigStimPlanFractureDefinition>>
RimFractureGroupStatistics::readFractureDefinitions( const std::vector<caf::FilePath>& filePaths,
RiaDefines::EclipseUnitSystem unitSystem ) const
{
double conductivityScaleFactor = 1.0;
std::vector<cvf::ref<RigStimPlanFractureDefinition>> results;
for ( auto filePath : m_filePaths.v() )
{
RiaLogging::info( QString( "Loading file: %1" ).arg( filePath.path() ) );
QString errorMessage;
cvf::ref<RigStimPlanFractureDefinition> stimPlanFractureDefinitionData =
RifStimPlanXmlReader::readStimPlanXMLFile( filePath.path(),
conductivityScaleFactor,
RifStimPlanXmlReader::MirrorMode::MIRROR_AUTO,
unitSystem,
&errorMessage );
if ( !errorMessage.isEmpty() )
{
RiaLogging::error( QString( "Error when reading file: '%1'" ).arg( errorMessage ) );
}
if ( stimPlanFractureDefinitionData.notNull() )
{
results.push_back( stimPlanFractureDefinitionData );
}
}
return results;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::set<QString> RimFractureGroupStatistics::findAllResultNames(
const std::vector<cvf::ref<RigStimPlanFractureDefinition>>& stimPlanFractureDefinitions )
{
std::set<QString> resultNames;
for ( auto stimPlanFractureDefinitionData : stimPlanFractureDefinitions )
{
for ( auto resultName : stimPlanFractureDefinitionData->conductivityResultNames() )
{
resultNames.insert( resultName );
}
}
return resultNames;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::vector<cvf::cref<RigFractureGrid>> RimFractureGroupStatistics::createFractureGrids(
const std::vector<cvf::ref<RigStimPlanFractureDefinition>>& stimPlanFractureDefinitions,
RiaDefines::EclipseUnitSystem unitSystem )
{
// Defaults to avoid scaling
double halfLengthScaleFactor = 1.0;
double heightScaleFactor = 1.0;
std::vector<cvf::cref<RigFractureGrid>> fractureGrids;
for ( auto stimPlanFractureDefinitionData : stimPlanFractureDefinitions )
{
QStringList conductivityResultNames = stimPlanFractureDefinitionData->conductivityResultNames();
double wellPathDepthAtFracture = computeDepthOfWellPathAtFracture( stimPlanFractureDefinitionData );
// Always use last time steps
std::vector<double> timeSteps = stimPlanFractureDefinitionData->timeSteps();
int activeTimeStepIndex = timeSteps.size() - 1;
QString conductivityResultNameOnFile = conductivityResultNames.first();
RiaLogging::info( QString( "Using result: %1" ).arg( conductivityResultNameOnFile ) );
cvf::cref<RigFractureGrid> fractureGrid =
stimPlanFractureDefinitionData->createFractureGrid( conductivityResultNameOnFile,
activeTimeStepIndex,
halfLengthScaleFactor,
heightScaleFactor,
wellPathDepthAtFracture,
unitSystem );
if ( fractureGrid.notNull() )
{
fractureGrids.push_back( fractureGrid );
}
}
return fractureGrids;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::tuple<double, double, double, double>
RimFractureGroupStatistics::findExtentsOfGrids( const std::vector<cvf::cref<RigFractureGrid>>& fractureGrids )
{
// Find min and max extent of all the grids
double minX = std::numeric_limits<double>::max();
double maxX = -std::numeric_limits<double>::max();
double minY = std::numeric_limits<double>::max();
double maxY = -std::numeric_limits<double>::max();
for ( auto fractureGrid : fractureGrids )
{
for ( auto fractureCell : fractureGrid->fractureCells() )
{
for ( auto polygon : fractureCell.getPolygon() )
{
minX = std::min( minX, polygon.x() );
maxX = std::max( maxX, polygon.x() );
minY = std::min( minY, polygon.y() );
maxY = std::max( maxY, polygon.y() );
}
}
}
return std::make_tuple( minX, maxX, minY, maxY );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
bool RimFractureGroupStatistics::isCoordinateInsideFractureCell( double x, double y, const RigFractureCell& cell )
{
const cvf::Vec3d& minPoint = cell.getPolygon()[0];
const cvf::Vec3d& maxPoint = cell.getPolygon()[2];
// TODO: Investigate strange ordering for y coords.
return ( x > minPoint.x() && x <= maxPoint.x() && y <= minPoint.y() && y > maxPoint.y() );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
double RimFractureGroupStatistics::computeDepthOfWellPathAtFracture(
cvf::ref<RigStimPlanFractureDefinition> stimPlanFractureDefinitionData )
{
double firstTvd = stimPlanFractureDefinitionData->topPerfTvd();
double lastTvd = stimPlanFractureDefinitionData->bottomPerfTvd();
if ( firstTvd != HUGE_VAL && lastTvd != HUGE_VAL )
{
return ( firstTvd + lastTvd ) / 2;
}
else
{
firstTvd = stimPlanFractureDefinitionData->minDepth();
lastTvd = stimPlanFractureDefinitionData->maxDepth();
return ( firstTvd + lastTvd ) / 2;
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RimFractureGroupStatistics::sampleAllGrids( const std::vector<cvf::cref<RigFractureGrid>>& fractureGrids,
std::vector<std::vector<double>>& samples,
double minX,
double minY,
int numSamplesX,
int numSamplesY,
double sampleDistanceX,
double sampleDistanceY )
{
for ( int y = 0; y < numSamplesY; y++ )
{
for ( int x = 0; x < numSamplesX; x++ )
{
double posX = minX + x * sampleDistanceX;
double posY = minY + y * sampleDistanceY;
for ( auto fractureGrid : fractureGrids )
{
for ( auto fractureCell : fractureGrid->fractureCells() )
{
if ( isCoordinateInsideFractureCell( posX, posY, fractureCell ) )
{
int idx = y * numSamplesX + x;
samples[idx].push_back( fractureCell.getConductivityValue() );
break;
}
}
}
}
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
bool RimFractureGroupStatistics::writeStatisticsToCsv( const QString& filePath, const RigSlice2D& samples )
{
QFile data( filePath );
if ( !data.open( QFile::WriteOnly | QFile::Truncate ) )
{
return false;
}
QTextStream stream( &data );
QString fieldSeparator = RiaPreferences::current()->csvTextExportFieldSeparator;
RifCsvDataTableFormatter formatter( stream, fieldSeparator );
std::vector<RifTextDataTableColumn> header;
for ( size_t y = 0; y < samples.ny(); y++ )
header.push_back( RifTextDataTableColumn( "", RifTextDataTableDoubleFormat::RIF_FLOAT ) );
formatter.header( header );
for ( size_t y = 0; y < samples.ny(); y++ )
{
for ( size_t x = 0; x < samples.nx(); x++ )
{
formatter.add( samples.getValue( x, y ) );
}
formatter.rowCompleted();
}
formatter.tableCompleted();
return true;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RimFractureGroupStatistics::generateStatisticsGrids( const std::vector<std::vector<double>>& samples,
int numSamplesX,
int numSamplesY,
std::map<QString, std::shared_ptr<RigSlice2D>>& statisticsGrids )
{
std::vector<QString> names = { "mean", "max", "min", "p10", "p50", "p90", "occurrence" };
for ( auto name : names )
{
std::shared_ptr<RigSlice2D> grid = std::make_shared<RigSlice2D>( numSamplesX, numSamplesY );
statisticsGrids[name] = grid;
}
for ( int y = 0; y < numSamplesY; y++ )
{
for ( int x = 0; x < numSamplesX; x++ )
{
size_t idx = y * numSamplesX + x;
double min;
double max;
double sum;
double range;
double mean;
double dev;
RigStatisticsMath::calculateBasicStatistics( samples[idx], &min, &max, &sum, &range, &mean, &dev );
statisticsGrids["mean"]->setValue( x, y, mean );
statisticsGrids["min"]->setValue( x, y, min );
statisticsGrids["max"]->setValue( x, y, max );
double p10;
double p50;
double p90;
RigStatisticsMath::calculateStatisticsCurves( samples[idx], &p10, &p50, &p90, &mean );
statisticsGrids["p10"]->setValue( x, y, p10 );
statisticsGrids["p50"]->setValue( x, y, p50 );
statisticsGrids["p90"]->setValue( x, y, p90 );
statisticsGrids["occurrence"]->setValue( x, y, samples[idx].size() );
}
}
}

View File

@ -20,6 +20,11 @@
#include "RimNamedObject.h"
#include "RigStimPlanFractureDefinition.h"
class RigFractureCell;
class RigSlice2D;
//==================================================================================================
///
///
@ -40,6 +45,40 @@ protected:
caf::PdmUiEditorAttribute* attribute ) override;
QString generateFilePathsTable();
void computeStatistics();
static double computeDepthOfWellPathAtFracture( cvf::ref<RigStimPlanFractureDefinition> stimPlanFractureDefinitionData );
static bool isCoordinateInsideFractureCell( double x, double y, const RigFractureCell& cell );
std::vector<cvf::ref<RigStimPlanFractureDefinition>>
readFractureDefinitions( const std::vector<caf::FilePath>& filePaths,
RiaDefines::EclipseUnitSystem unitSystem ) const;
std::vector<cvf::cref<RigFractureGrid>>
createFractureGrids( const std::vector<cvf::ref<RigStimPlanFractureDefinition>>& stimPlanFractureDefinitions,
RiaDefines::EclipseUnitSystem unitSystem );
static std::set<QString>
findAllResultNames( const std::vector<cvf::ref<RigStimPlanFractureDefinition>>& stimPlanFractureDefinitions );
static std::tuple<double, double, double, double>
findExtentsOfGrids( const std::vector<cvf::cref<RigFractureGrid>>& fractureGrids );
static void sampleAllGrids( const std::vector<cvf::cref<RigFractureGrid>>& fractureGrids,
std::vector<std::vector<double>>& samples,
double minX,
double maxX,
int numSamplesX,
int numSamplesY,
double sampleDistanceX,
double sampleDistanceY );
static void generateStatisticsGrids( const std::vector<std::vector<double>>& samples,
int numSamplesX,
int numSamplesY,
std::map<QString, std::shared_ptr<RigSlice2D>>& statisticsGrids );
static bool writeStatisticsToCsv( const QString& filePath, const RigSlice2D& samples );
caf::PdmField<std::vector<caf::FilePath>> m_filePaths;
caf::PdmField<QString> m_filePathsTable;
};