#3201 #3202 Split in separate files. Move offshore spherical coords stuff to separate file in RiaTools

This commit is contained in:
Jacob Støren
2018-08-13 16:17:31 +02:00
parent c6069ade29
commit fb7a7e39eb
12 changed files with 390 additions and 113 deletions

View File

@@ -30,6 +30,9 @@ ${CMAKE_CURRENT_LIST_DIR}/RiaTimeHistoryCurveResampler.h
${CMAKE_CURRENT_LIST_DIR}/RiaStatisticsTools.h
${CMAKE_CURRENT_LIST_DIR}/RiaPolyArcLineSampler.h
${CMAKE_CURRENT_LIST_DIR}/RiaSCurveCalculator.h
${CMAKE_CURRENT_LIST_DIR}/RiaArcCurveCalculator.h
${CMAKE_CURRENT_LIST_DIR}/RiaJCurveCalculator.h
${CMAKE_CURRENT_LIST_DIR}/RiaOffshoreSphericalCoords.h
)
set (SOURCE_GROUP_SOURCE_FILES
@@ -63,6 +66,8 @@ ${CMAKE_CURRENT_LIST_DIR}/RiaTimeHistoryCurveResampler.cpp
${CMAKE_CURRENT_LIST_DIR}/RiaStatisticsTools.cpp
${CMAKE_CURRENT_LIST_DIR}/RiaPolyArcLineSampler.cpp
${CMAKE_CURRENT_LIST_DIR}/RiaSCurveCalculator.cpp
${CMAKE_CURRENT_LIST_DIR}/RiaArcCurveCalculator.cpp
${CMAKE_CURRENT_LIST_DIR}/RiaJCurveCalculator.cpp
)
list(APPEND CODE_HEADER_FILES

View File

@@ -0,0 +1,79 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2018- 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 "RiaArcCurveCalculator.h"
#include "RiaOffshoreSphericalCoords.h"
//--------------------------------------------------------------------------------------------------
/// + p1
/// t1 //
/// | + C
/// \
/// + p2
//--------------------------------------------------------------------------------------------------
RiaArcCurveCalculator::RiaArcCurveCalculator(cvf::Vec3d p1, cvf::Vec3d t1, cvf::Vec3d p2)
: m_isCalculationOK(false)
, m_radius(std::numeric_limits<double>::infinity())
, m_arcCS(cvf::Mat4d::ZERO)
{
bool isOk = t1.normalize();
if (!isOk)
{
// No tangent. Bail out
return;
}
cvf::Vec3d p1p2 = p2 - p1;
cvf::Vec3d t12 = p1p2.getNormalized(&isOk);
if (!isOk)
{
// p1 and p2 in the same place.
return;
}
cvf::Vec3d N = (t1 ^ t12).getNormalized(&isOk);
if (!isOk)
{
// P2 is on the p1 + k*t1 line. We have a straight line
return;
}
cvf::Vec3d tr1 = (N ^ t1).getNormalized();
m_radius = 0.5 * p1p2.length() / (tr1.dot(t12));
cvf::Vec3d C = p1 + m_radius * tr1;
cvf::Vec3d nTr1 = -tr1;
m_arcCS = cvf::Mat4d::fromCoordSystemAxes( &nTr1, &t1, &N );
m_arcCS.setTranslation(C);
m_isCalculationOK = true;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RiaArcCurveCalculator::RiaArcCurveCalculator(cvf::Vec3d p1, double azi1, double inc1, cvf::Vec3d p2)
{
cvf::Vec3d t1( RiaOffshoreSphericalCoords::unitVectorFromAziInc(azi1,inc1));
(*this) = RiaArcCurveCalculator(p1, t1, p2);
}

View File

@@ -0,0 +1,51 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2018- 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 "cvfBase.h"
#include "cvfVector3.h"
#include "cvfMatrix4.h"
//--------------------------------------------------------------------------------------------------
/// + p1
/// t1 //
/// | + C
/// \
/// + p2
//--------------------------------------------------------------------------------------------------
class RiaArcCurveCalculator
{
public:
RiaArcCurveCalculator(cvf::Vec3d p1, cvf::Vec3d t1, cvf::Vec3d p2);
RiaArcCurveCalculator(cvf::Vec3d p1, double azi1, double inc1, cvf::Vec3d p2);
bool isOk() const { return m_isCalculationOK;}
cvf::Mat4d arcCS() const { return m_arcCS; }
double radius() const { return m_radius;}
cvf::Vec3d center() const { return m_arcCS.translation();}
cvf::Vec3d normal() const { return cvf::Vec3d(m_arcCS.col(2));}
private:
bool m_isCalculationOK;
double m_radius;
cvf::Mat4d m_arcCS;
};

View File

@@ -0,0 +1,70 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2018- 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 "RiaJCurveCalculator.h"
#include "RiaOffshoreSphericalCoords.h"
#include "cvfMatrix3.h"
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RiaJCurveCalculator::RiaJCurveCalculator(cvf::Vec3d p1, double azi1, double inc1, double r1,
cvf::Vec3d p2)
: m_isCalculationOK(false)
{
cvf::Vec3d t1 (RiaOffshoreSphericalCoords::unitVectorFromAziInc(azi1, inc1));
cvf::Vec3d p1p2 = p2 - p1;
bool isOk = true;
cvf::Vec3d tr1 = (p1p2 - (p1p2.dot(t1)) * t1).getNormalized(&isOk);
if (!isOk)
{
// p2 is on the p1 + t12 line. Degenerates to a line.
m_firstArcEndpoint = p2;
m_c1 = cvf::Vec3d::UNDEFINED;
m_n1 = cvf::Vec3d::UNDEFINED;
return;
}
cvf::Vec3d c1 = p1 + r1 * tr1;
cvf::Vec3d p2c1 = c1 - p2;
double p2c1Length = p2c1.length();
if (p2c1Length < r1)
{
// Radius is too big. We can not get to point 2 using the requested radius.
}
double d = sqrt( p2c1Length * p2c1Length - r1 * r1);
double betha = asin( r1/p2c1Length );
cvf::Vec3d tp2c1 = p2c1/p2c1Length;
cvf::Vec3d nc1 = t1 ^ tr1;
cvf::Vec3d tp11p2 = tp2c1.getTransformedVector(cvf::Mat3d::fromRotation(nc1, betha));
m_firstArcEndpoint = p2 - d*tp11p2;
m_c1 = c1;
m_n1 = nc1;
m_isCalculationOK = true;
}

View File

@@ -0,0 +1,43 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2018- 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 "cvfBase.h"
#include "cvfVector3.h"
class RiaJCurveCalculator
{
public:
RiaJCurveCalculator(cvf::Vec3d p1, double azi1, double inc1, double r1,
cvf::Vec3d p2);
bool isOk() const { return m_isCalculationOK;}
cvf::Vec3d firstArcEndpoint() const { return m_firstArcEndpoint; }
cvf::Vec3d firstCenter() const { return m_c1; }
cvf::Vec3d firstNormal() const { return m_n1; }
private:
bool m_isCalculationOK;
cvf::Vec3d m_firstArcEndpoint;
cvf::Vec3d m_c1;
cvf::Vec3d m_n1;
};

View File

@@ -0,0 +1,62 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2018- 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 <array>
#include "cvfVector3.h"
#include <cmath>
// Y - North, X - East, Z - up but depth is negative Z
// azi is measured from the Northing (Y) Axis in Clockwise direction looking down
// inc is measured from the negative Z (depth) axis
class RiaOffshoreSphericalCoords
{
public:
explicit RiaOffshoreSphericalCoords(const cvf::Vec3f& vec)
{
// Azimuth:
if (vec[0] == 0.0f && vec[1] == 0.0 ) incAziR[1] = 0.0f;
else incAziR[1] = atan2(vec[0], vec[1]); // atan2(Y, X)
// R
incAziR[2] = vec.length();
// Inclination from vertical down
if (incAziR[2] == 0) incAziR[0] = 0.0f;
else incAziR[0] = acos(-vec[2]/incAziR[2]);
}
float inc() const { return incAziR[0];}
float azi() const { return incAziR[1];}
float r() const { return incAziR[2];}
// Note that this is a double function, while the rest of the class is float.
// Todo: Convert class to a template to enable float and double versions of everything
static cvf::Vec3d unitVectorFromAziInc(double azimuth, double inclination)
{
return cvf::Vec3d(sin(azimuth)*sin(inclination),
cos(azimuth)*sin(inclination),
-cos(inclination));
}
private:
std::array<float, 3> incAziR;
};

View File

@@ -1,5 +1,25 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2018- 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 "RiaSCurveCalculator.h"
#include "RiaOffshoreSphericalCoords.h"
#include "SolveSpaceSystem.h"
#include <cmath>
#include "cvfMatrix4.h"
@@ -78,12 +98,8 @@ RiaSCurveCalculator::RiaSCurveCalculator(cvf::Vec3d p1, double azi1, double inc1
estimatedCurveCalc.dump();
#endif
cvf::Vec3d t1(sin(azi1)*sin(inc1),
cos(azi1)*sin(inc1),
-cos(inc1));
cvf::Vec3d t2(sin(azi2)*sin(inc2),
cos(azi2)*sin(inc2),
-cos(inc2));
cvf::Vec3d t1(RiaOffshoreSphericalCoords::unitVectorFromAziInc(azi1,inc1));
cvf::Vec3d t2(RiaOffshoreSphericalCoords::unitVectorFromAziInc(azi2,inc2));
cvf::Vec3d est_tp1c1 = (estimatedCurveCalc.firstCenter() - p1).getNormalized();
cvf::Vec3d est_tp2c2 = (estimatedCurveCalc.secondCenter() - p2).getNormalized();
@@ -638,7 +654,7 @@ RiaSCurveCalculator::RiaSCurveCalculator(cvf::Vec3d p1, cvf::Vec3d q1,
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RiaSCurveCalculator::dump()
void RiaSCurveCalculator::dump() const
{
cvf::Vec3d v_C1 = firstCenter();
cvf::Vec3d v_C2 = secondCenter();
@@ -666,12 +682,8 @@ void RiaSCurveCalculator::dump()
RiaSCurveCalculator RiaSCurveCalculator::fromTangentsAndLength(cvf::Vec3d p1, double azi1, double inc1, double lengthToQ1,
cvf::Vec3d p2, double azi2, double inc2, double lengthToQ2)
{
cvf::Vec3d t1(sin(azi1)*sin(inc1),
cos(azi1)*sin(inc1),
-cos(inc1));
cvf::Vec3d t2(sin(azi2)*sin(inc2),
cos(azi2)*sin(inc2),
-cos(inc2));
cvf::Vec3d t1(RiaOffshoreSphericalCoords::unitVectorFromAziInc(azi1,inc1));
cvf::Vec3d t2(RiaOffshoreSphericalCoords::unitVectorFromAziInc(azi2,inc2));
cvf::Vec3d Q1 = p1 + lengthToQ1 * t1;
cvf::Vec3d Q2 = p2 - lengthToQ2 * t2;
@@ -682,16 +694,3 @@ RiaSCurveCalculator RiaSCurveCalculator::fromTangentsAndLength(cvf::Vec3d p1, do
return curveFromControlPoints;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RiaSCurveCalculator::calculateEstimatedSolution()
{
// Plane1 basisvectors
// C1 position in Plane1
// P11 position in Plane 1
// Plane2 basisvectors
}

View File

@@ -1,4 +1,22 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2018- 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 "cvfBase.h"
#include "cvfVector3.h"
@@ -11,52 +29,32 @@ public:
RiaSCurveCalculator( cvf::Vec3d p1, cvf::Vec3d q1,
cvf::Vec3d p2, cvf::Vec3d q2 );
bool isOk() { return m_isCalculationOK;}
cvf::Vec3d firstArcEndpoint() { return m_firstArcEndpoint; }
cvf::Vec3d secondArcStartpoint() { return m_secondArcStartpoint; }
cvf::Vec3d firstCenter() { return m_c1; }
cvf::Vec3d secondCenter() { return m_c2; }
cvf::Vec3d firstNormal() { return m_n1; }
cvf::Vec3d secondNormal() { return m_n2; }
double firstRadius() { return m_r1; }
double secondRadius() { return m_r2; }
bool isOk() const { return m_isCalculationOK; }
cvf::Vec3d firstArcEndpoint() const { return m_firstArcEndpoint; }
cvf::Vec3d secondArcStartpoint() const { return m_secondArcStartpoint; }
cvf::Vec3d firstCenter() const { return m_c1; }
cvf::Vec3d secondCenter() const { return m_c2; }
cvf::Vec3d firstNormal() const { return m_n1; }
cvf::Vec3d secondNormal() const { return m_n2; }
double firstRadius() const { return m_r1; }
double secondRadius() const { return m_r2; }
void dump();
void dump() const;
static RiaSCurveCalculator fromTangentsAndLength(cvf::Vec3d p1, double azi1, double inc1, double lengthToQ1,
cvf::Vec3d p2, double azi2, double inc2, double lengthToQ2 );
private:
void calculateEstimatedSolution();
private:
bool m_isCalculationOK;
cvf::Vec3d m_p1;
cvf::Vec3d m_p2;
cvf::Vec3d m_firstArcEndpoint;
cvf::Vec3d m_secondArcStartpoint;
cvf::Vec3d m_c1;
cvf::Vec3d m_c2;
cvf::Vec3d m_n1;
cvf::Vec3d m_n2;
double m_r1;
double m_r2;
double m_r1;
double m_r2;
};