Faults: Reading from Eclipse text files and basic visualization

Added reading of faults directly from Eclipse text files
Added named faults in project tree
Added fault geompetry parts
Added test data
This commit is contained in:
Magne Sjaastad
2013-12-03 20:30:32 +01:00
parent a9140ca894
commit 12918ec107
39 changed files with 2468 additions and 109 deletions

View File

@@ -0,0 +1,304 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2011-2012 Statoil ASA, Ceetron 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 "RivFaultGeometryGenerator.h"
#include "cvfDrawableGeo.h"
#include "cvfPrimitiveSetIndexedUInt.h"
#include "cvfOutlineEdgeExtractor.h"
#include "cvfStructGridScalarDataAccess.h"
#include "cvfScalarMapper.h"
#include "RigFault.h"
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RivFaultGeometryGenerator::RivFaultGeometryGenerator(const cvf::StructGridInterface* grid, const RigFault* fault)
: m_grid(grid),
m_fault(fault)
{
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RivFaultGeometryGenerator::~RivFaultGeometryGenerator()
{
}
//--------------------------------------------------------------------------------------------------
/// Generate surface drawable geo from the specified region
//--------------------------------------------------------------------------------------------------
cvf::ref<cvf::DrawableGeo> RivFaultGeometryGenerator::generateSurface()
{
computeArrays();
CVF_ASSERT(m_vertices.notNull());
if (m_vertices->size() == 0) return NULL;
cvf::ref<cvf::DrawableGeo> geo = new cvf::DrawableGeo;
geo->setFromQuadVertexArray(m_vertices.p());
return geo;
}
//--------------------------------------------------------------------------------------------------
/// Generates simplified mesh as line drawing
/// Must call generateSurface first
//--------------------------------------------------------------------------------------------------
cvf::ref<cvf::DrawableGeo> RivFaultGeometryGenerator::createMeshDrawable()
{
if (!(m_vertices.notNull() && m_vertices->size() != 0)) return NULL;
cvf::ref<cvf::DrawableGeo> geo = new cvf::DrawableGeo;
geo->setVertexArray(m_vertices.p());
cvf::ref<cvf::UIntArray> indices = lineIndicesFromQuadVertexArray(m_vertices.p());
cvf::ref<cvf::PrimitiveSetIndexedUInt> prim = new cvf::PrimitiveSetIndexedUInt(cvf::PT_LINES);
prim->setIndices(indices.p());
geo->addPrimitiveSet(prim.p());
return geo;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
cvf::ref<cvf::DrawableGeo> RivFaultGeometryGenerator::createOutlineMeshDrawable(double creaseAngle)
{
if (!(m_vertices.notNull() && m_vertices->size() != 0)) return NULL;
cvf::OutlineEdgeExtractor ee(creaseAngle, *m_vertices);
cvf::ref<cvf::UIntArray> indices = lineIndicesFromQuadVertexArray(m_vertices.p());
ee.addPrimitives(4, *indices);
cvf::ref<cvf::UIntArray> lineIndices = ee.lineIndices();
if (lineIndices->size() == 0)
{
return NULL;
}
cvf::ref<cvf::PrimitiveSetIndexedUInt> prim = new cvf::PrimitiveSetIndexedUInt(cvf::PT_LINES);
prim->setIndices(lineIndices.p());
cvf::ref<cvf::DrawableGeo> geo = new cvf::DrawableGeo;
geo->setVertexArray(m_vertices.p());
geo->addPrimitiveSet(prim.p());
return geo;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
cvf::ref<cvf::UIntArray> RivFaultGeometryGenerator::lineIndicesFromQuadVertexArray(const cvf::Vec3fArray* vertexArray)
{
CVF_ASSERT(vertexArray);
size_t numVertices = vertexArray->size();
int numQuads = static_cast<int>(numVertices/4);
CVF_ASSERT(numVertices%4 == 0);
cvf::ref<cvf::UIntArray> indices = new cvf::UIntArray;
indices->resize(numQuads*8);
#pragma omp parallel for
for (int i = 0; i < numQuads; i++)
{
int idx = 8*i;
indices->set(idx + 0, i*4 + 0);
indices->set(idx + 1, i*4 + 1);
indices->set(idx + 2, i*4 + 1);
indices->set(idx + 3, i*4 + 2);
indices->set(idx + 4, i*4 + 2);
indices->set(idx + 5, i*4 + 3);
indices->set(idx + 6, i*4 + 3);
indices->set(idx + 7, i*4 + 0);
}
return indices;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RivFaultGeometryGenerator::computeArrays()
{
std::vector<cvf::Vec3f> vertices;
m_quadsToGridCells.clear();
m_quadsToFace.clear();
cvf::Vec3d offset = m_grid->displayModelOffset();
for (size_t faceType = 0; faceType < 6; faceType++)
{
cvf::StructGridInterface::FaceType faceEnum = cvf::StructGridInterface::FaceType(faceType);
const std::vector<cvf::CellRange>& cellRanges = m_fault->cellRangeForFace(faceEnum);
for (size_t i = 0; i < cellRanges.size(); i++)
{
const cvf::CellRange& cellRange = cellRanges[i];
std::vector<size_t> gridCellIndices;
computeGlobalCellIndices(cellRange, gridCellIndices);
for (size_t gIdx = 0; gIdx < gridCellIndices.size(); gIdx++)
{
size_t cellIndex = gridCellIndices[gIdx];
if (!(*m_cellVisibility)[cellIndex]) continue;
cvf::Vec3d cornerVerts[8];
m_grid->cellCornerVertices(cellIndex, cornerVerts);
cvf::ubyte faceConn[4];
m_grid->cellFaceVertexIndices(faceEnum, faceConn);
// Critical section to avoid two threads accessing the arrays at the same time.
#pragma omp critical
{
int n;
for (n = 0; n < 4; n++)
{
vertices.push_back(cvf::Vec3f(cornerVerts[faceConn[n]] - offset));
}
// Keep track of the source cell index per quad
m_quadsToGridCells.push_back(cellIndex);
m_quadsToFace.push_back(faceEnum);
}
}
}
}
m_vertices = new cvf::Vec3fArray;
m_vertices->assign(vertices);
}
//--------------------------------------------------------------------------------------------------
/// Calculates the texture coordinates in a "nearly" one dimensional texture.
/// Undefined values are coded with a y-texture coordinate value of 1.0 instead of the normal 0.5
//--------------------------------------------------------------------------------------------------
void RivFaultGeometryGenerator::textureCoordinates(cvf::Vec2fArray* textureCoords, const cvf::StructGridScalarDataAccess* dataAccessObject, const cvf::ScalarMapper* mapper) const
{
if (!dataAccessObject) return;
size_t numVertices = m_quadsToGridCells.size()*4;
textureCoords->resize(numVertices);
cvf::Vec2f* rawPtr = textureCoords->ptr();
double cellScalarValue;
cvf::Vec2f texCoord;
#pragma omp parallel for private(texCoord, cellScalarValue)
for (int i = 0; i < static_cast<int>(m_quadsToGridCells.size()); i++)
{
cellScalarValue = dataAccessObject->cellScalar(m_quadsToGridCells[i]);
texCoord = mapper->mapToTextureCoord(cellScalarValue);
if (cellScalarValue == HUGE_VAL || cellScalarValue != cellScalarValue) // a != a is true for NAN's
{
texCoord[1] = 1.0f;
}
size_t j;
for (j = 0; j < 4; j++)
{
rawPtr[i*4 + j] = texCoord;
}
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
cvf::ref<cvf::Array<size_t> > RivFaultGeometryGenerator::triangleToSourceGridCellMap() const
{
cvf::ref<cvf::Array<size_t> > triangles = new cvf::Array<size_t>(2*m_quadsToGridCells.size());
#pragma omp parallel for
for (int i = 0; i < static_cast<int>(m_quadsToGridCells.size()); i++)
{
triangles->set(i*2, m_quadsToGridCells[i]);
triangles->set(i*2+1, m_quadsToGridCells[i]);
}
return triangles;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RivFaultGeometryGenerator::setCellVisibility(const cvf::UByteArray* cellVisibility)
{
m_cellVisibility = cellVisibility;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<size_t>& RivFaultGeometryGenerator::quadToGridCellIndices() const
{
return m_quadsToGridCells;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<cvf::StructGridInterface::FaceType>& RivFaultGeometryGenerator::quadToFace() const
{
return m_quadsToFace;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RivFaultGeometryGenerator::computeGlobalCellIndices(const cvf::CellRange& cellRange, std::vector<size_t>& gridCellIndices) const
{
cvf::Vec3st min, max;
cellRange.range(min, max);
for (size_t i = min.x(); i <= max.x(); i++)
{
if (i >= m_grid->cellCountI())
{
continue;
}
for (size_t j = min.y(); j <= max.y(); j++)
{
if (j >= m_grid->cellCountJ())
{
continue;
}
for (size_t k = min.z(); k <= max.z(); k++)
{
if (k >= m_grid->cellCountK())
{
continue;
}
gridCellIndices.push_back(m_grid->cellIndexFromIJK(i, j, k));
}
}
}
}