2 #include "HYDROData_Bathymetry.h"
3 #include "HYDROData_Document.h"
4 #include "HYDROData_Tool.h"
9 #include <TDataStd_RealArray.hxx>
10 #include <TDataStd_AsciiString.hxx>
16 #include <QStringList>
18 #define INVALID_ALTITUDE_VALUE -9999.0
19 #define PYTHON_BATHYMETRY_ID "KIND_BATHYMETRY"
22 IMPLEMENT_STANDARD_HANDLE(HYDROData_Bathymetry, HYDROData_Entity)
23 IMPLEMENT_STANDARD_RTTIEXT(HYDROData_Bathymetry, HYDROData_Entity)
25 HYDROData_Bathymetry::HYDROData_Bathymetry()
30 HYDROData_Bathymetry::~HYDROData_Bathymetry()
34 QStringList HYDROData_Bathymetry::DumpToPython( MapOfTreatedObjects& theTreatedObjects ) const
38 Handle(HYDROData_Document) aDocument = HYDROData_Document::Document( this );
39 if ( aDocument.IsNull() )
42 QString aDocName = aDocument->GetDocPyName();
43 QString aBathymetryName = GetName();
45 aResList << QString( "%1 = %2.CreateObject( %3 );" )
46 .arg( aBathymetryName ).arg( aDocName ).arg( PYTHON_BATHYMETRY_ID );
47 aResList << QString( "%1.SetName( \"%2\" );" )
48 .arg( aBathymetryName ).arg( aBathymetryName );
50 QString aFilePath = GetFilePath();
51 if ( !aFilePath.isEmpty() )
53 aResList << QString( "%1.ImportFromFile( \"%2\" );" )
54 .arg( aBathymetryName ).arg( aFilePath );
58 // TODO : bathymetry is composed from other bathymetry(ies)
64 double HYDROData_Bathymetry::GetInvalidAltitude()
66 return INVALID_ALTITUDE_VALUE;
69 void HYDROData_Bathymetry::SetAltitudePoints( const AltitudePoints& thePoints )
71 RemoveAltitudePoints();
73 if ( thePoints.isEmpty() )
77 Handle(TDataStd_RealArray) aCoordsArray =
78 TDataStd_RealArray::Set( myLab.FindChild( DataTag_AltitudePoints ), 0, thePoints.size() * 3 - 1 );
80 AltitudePoints::const_iterator aListItBeg = thePoints.constBegin();
81 AltitudePoints::const_iterator aListItEnd = thePoints.constEnd();
82 for ( int i = 0 ; aListItBeg != aListItEnd; ++i, ++aListItBeg )
84 const AltitudePoint& aPoint = *aListItBeg;
86 aCoordsArray->SetValue( i * 3, aPoint.X() );
87 aCoordsArray->SetValue( i * 3 + 1, aPoint.Y() );
88 aCoordsArray->SetValue( i * 3 + 2, aPoint.Z() );
92 HYDROData_Bathymetry::AltitudePoints HYDROData_Bathymetry::GetAltitudePoints() const
94 AltitudePoints aPoints;
96 Handle(TDataStd_RealArray) aCoordsArray;
97 if ( !myLab.FindChild( DataTag_AltitudePoints ).FindAttribute( TDataStd_RealArray::GetID(), aCoordsArray ) )
100 int aLowerIdx = aCoordsArray->Lower();
101 int anUpperIdx = aCoordsArray->Upper();
102 for ( int i = aCoordsArray->Lower(), n = aCoordsArray->Upper(); i <= n; )
107 AltitudePoint aPoint;
108 aPoint.SetX( aCoordsArray->Value( i++ ) );
109 aPoint.SetY( aCoordsArray->Value( i++ ) );
110 aPoint.SetZ( aCoordsArray->Value( i++ ) );
118 void HYDROData_Bathymetry::RemoveAltitudePoints()
120 TDF_Label aLab = myLab.FindChild( DataTag_AltitudePoints );
121 aLab.ForgetAllAttributes();
124 void interpolateAltitudeForPoints( const gp_XY& thePoint,
125 const HYDROData_Bathymetry::AltitudePoint& theFirstPoint,
126 const HYDROData_Bathymetry::AltitudePoint& theSecPoint,
127 HYDROData_Bathymetry::AltitudePoint& theResPoint,
128 const bool& theIsVertical )
130 double aCoordX = thePoint.X();
131 double aCoordY = thePoint.Y();
135 aCoordX = theFirstPoint.X();
137 if ( !ValuesEquals( theFirstPoint.X(), theSecPoint.X() ) )
139 // Recalculate X coordinate by equation of line from two points
140 aCoordX = ( ( ( thePoint.Y() - theFirstPoint.Y() ) * ( theSecPoint.X() - theFirstPoint.X() ) ) /
141 ( theSecPoint.Y() - theFirstPoint.Y() ) ) + theFirstPoint.X();
146 aCoordY = theFirstPoint.Y();
148 if ( !ValuesEquals( theFirstPoint.Y(), theSecPoint.Y() ) )
150 // Recalculate y by equation of line from two points
151 aCoordY = ( ( ( thePoint.X() - theFirstPoint.X() ) * ( theSecPoint.Y() - theFirstPoint.Y() ) ) /
152 ( theSecPoint.X() - theFirstPoint.X() ) ) + theFirstPoint.Y();
156 theResPoint.SetX( aCoordX );
157 theResPoint.SetY( aCoordY );
159 // Calculate coefficient for interpolation
160 double aLength = Sqrt( Pow( theSecPoint.Y() - theFirstPoint.Y(), 2 ) +
161 Pow( theSecPoint.X() - theFirstPoint.X(), 2 ) );
163 double aInterCoeff = 0;
165 aInterCoeff = ( theSecPoint.Z() - theFirstPoint.Z() ) / aLength;
168 double aNewLength = Sqrt( Pow( theResPoint.Y() - theFirstPoint.Y(), 2 ) +
169 Pow( theResPoint.X() - theFirstPoint.X(), 2 ) );
171 // Calculate interpolated value
172 double aResVal = theFirstPoint.Z() + aInterCoeff * aNewLength;
174 theResPoint.SetZ( aResVal );
177 double HYDROData_Bathymetry::GetAltitudeForPoint( const QPointF& thePoint ) const
179 gp_XY aGpPoint( thePoint.x(), thePoint.y() );
180 return GetAltitudeForPoint( aGpPoint );
183 double HYDROData_Bathymetry::GetAltitudeForPoint( const gp_XY& thePoint ) const
185 double aResAltitude = -9999.90;
187 AltitudePoints anAltitudePoints = GetAltitudePoints();
188 if ( anAltitudePoints.isEmpty() )
191 QPolygonF aBoundingRect;
194 // [ 0 (top-left) ] [ 1 (top-right) ]
196 // [ 2 (bot-left) ] [ 3 (bot-right) ]
197 AltitudePoint aBounds[ 4 ] = { AltitudePoint( -DBL_MAX, -DBL_MAX, INVALID_ALTITUDE_VALUE ),
198 AltitudePoint( DBL_MAX, -DBL_MAX, INVALID_ALTITUDE_VALUE ),
199 AltitudePoint( -DBL_MAX, DBL_MAX, INVALID_ALTITUDE_VALUE ),
200 AltitudePoint( DBL_MAX, DBL_MAX, INVALID_ALTITUDE_VALUE ) };
202 AltitudePoints::const_iterator aListItBeg = anAltitudePoints.constBegin();
203 AltitudePoints::const_iterator aListItEnd = anAltitudePoints.constEnd();
204 for ( ; aListItBeg != aListItEnd; ++aListItBeg )
206 const AltitudePoint& aPoint = *aListItBeg;
208 double aDeltaX = Abs( aPoint.X() ) - Abs( thePoint.X() );
209 double aDeltaY = Abs( aPoint.Y() ) - Abs( thePoint.Y() );
211 if ( ValuesEquals( aDeltaX, 0.0 ) ) // Both left and right sides
213 if ( ValuesEquals( aDeltaY, 0.0 ) ) // Both top and bottom sides
215 aResAltitude = aPoint.Z();
218 else if ( aDeltaY < 0 ) // top side
221 if ( ValuesMoreEquals( aPoint.X(), aBounds[ 0 ].X() ) && ValuesMoreEquals( aPoint.Y(), aBounds[ 0 ].Y() ) )
222 aBounds[ 0 ] = aPoint;
223 if ( ValuesLessEquals( aPoint.X(), aBounds[ 1 ].X() ) && ValuesMoreEquals( aPoint.Y(), aBounds[ 1 ].Y() ) )
224 aBounds[ 1 ] = aPoint;
229 if ( ValuesMoreEquals( aPoint.X(), aBounds[ 2 ].X() ) && ValuesLessEquals( aPoint.Y(), aBounds[ 2 ].Y() ) )
230 aBounds[ 2 ] = aPoint;
231 if ( ValuesLessEquals( aPoint.X(), aBounds[ 3 ].X() ) && ValuesLessEquals( aPoint.Y(), aBounds[ 3 ].Y() ) )
232 aBounds[ 3 ] = aPoint;
235 else if ( aDeltaX < 0 ) // left side
237 if ( ValuesEquals( aDeltaY, 0.0 ) )
240 if ( ValuesMoreEquals( aPoint.X(), aBounds[ 0 ].X() ) && ValuesMoreEquals( aPoint.Y(), aBounds[ 0 ].Y() ) )
241 aBounds[ 0 ] = aPoint;
242 if ( ValuesMoreEquals( aPoint.X(), aBounds[ 2 ].X() ) && ValuesLessEquals( aPoint.Y(), aBounds[ 2 ].Y() ) )
243 aBounds[ 2 ] = aPoint;
245 else if ( aDeltaY < 0 )
248 if ( ValuesMoreEquals( aPoint.X(), aBounds[ 0 ].X() ) && ValuesMoreEquals( aPoint.Y(), aBounds[ 0 ].Y() ) )
249 aBounds[ 0 ] = aPoint;
253 // bottom left corner
254 if ( ValuesMoreEquals( aPoint.X(), aBounds[ 2 ].X() ) && ValuesLessEquals( aPoint.Y(), aBounds[ 2 ].Y() ) )
255 aBounds[ 2 ] = aPoint;
260 if ( ValuesEquals( aDeltaY, 0.0 ) )
263 if ( ValuesLessEquals( aPoint.X(), aBounds[ 1 ].X() ) && ValuesMoreEquals( aPoint.Y(), aBounds[ 1 ].Y() ) )
264 aBounds[ 1 ] = aPoint;
265 if ( ValuesLessEquals( aPoint.X(), aBounds[ 3 ].X() ) && ValuesLessEquals( aPoint.Y(), aBounds[ 3 ].Y() ) )
266 aBounds[ 3 ] = aPoint;
268 else if ( aDeltaY < 0 )
271 if ( ValuesLessEquals( aPoint.X(), aBounds[ 1 ].X() ) && ValuesMoreEquals( aPoint.Y(), aBounds[ 1 ].Y() ) )
272 aBounds[ 1 ] = aPoint;
276 // bottom right corner
277 if ( ValuesLessEquals( aPoint.X(), aBounds[ 3 ].X() ) && ValuesLessEquals( aPoint.Y(), aBounds[ 3 ].Y() ) )
278 aBounds[ 3 ] = aPoint;
282 // Update bounding rectangle of our global grid
283 aBoundingRect << QPointF( aPoint.X(), aPoint.Y() );
286 // Check if requested point is inside of our bounding rectangle
287 if ( !aBoundingRect.boundingRect().contains( thePoint.X(), thePoint.Y() ) )
288 return INVALID_ALTITUDE_VALUE;
290 // Calculate result altitude for point
291 AltitudePoint aFirstPoint( aBounds[ 0 ] ), aSecPoint( aBounds[ 1 ] );
293 // At first we merge top and bottom borders
294 if ( aBounds[ 0 ].Y() != aBounds[ 2 ].Y() || aBounds[ 0 ].X() != aBounds[ 2 ].X() )
295 interpolateAltitudeForPoints( thePoint, aBounds[ 0 ], aBounds[ 2 ], aFirstPoint, true );
297 if ( aBounds[ 1 ].Y() != aBounds[ 3 ].Y() || aBounds[ 1 ].X() != aBounds[ 3 ].X() )
298 interpolateAltitudeForPoints( thePoint, aBounds[ 1 ], aBounds[ 3 ], aSecPoint, true );
300 AltitudePoint aResPoint( aFirstPoint );
302 // At last we merge left and right borders
303 if ( aFirstPoint.Y() != aSecPoint.Y() || aFirstPoint.X() != aSecPoint.X() )
304 interpolateAltitudeForPoints( thePoint, aFirstPoint, aSecPoint, aResPoint, false );
306 aResAltitude = aResPoint.Z();
311 void HYDROData_Bathymetry::SetFilePath(const QString& theFilePath)
313 TCollection_AsciiString anAsciiStr( theFilePath.toStdString().c_str() );
314 TDataStd_AsciiString::Set( myLab.FindChild( DataTag_FilePath ), anAsciiStr );
317 QString HYDROData_Bathymetry::GetFilePath() const
321 Handle(TDataStd_AsciiString) anAsciiStr;
322 if ( myLab.FindChild( DataTag_FilePath ).FindAttribute( TDataStd_AsciiString::GetID(), anAsciiStr ) )
323 aRes = QString( anAsciiStr->Get().ToCString() );
328 bool HYDROData_Bathymetry::ImportFromFile( const QString& theFileName )
330 // Try to open the file
331 QFile aFile( theFileName );
332 if ( !aFile.exists() || !aFile.open( QIODevice::ReadOnly ) )
337 QString aFileSuf = QFileInfo( aFile ).suffix().toLower();
339 AltitudePoints aPoints;
341 // Try to import the file
342 if ( aFileSuf == "xyz" )
343 aRes = importFromXYZFile( aFile, aPoints );
350 // Update file path and altitude points of this Bathymetry
351 SetFilePath( theFileName );
352 SetAltitudePoints( aPoints );
355 return aRes && !aPoints.isEmpty();
358 bool HYDROData_Bathymetry::importFromXYZFile( QFile& theFile,
359 AltitudePoints& thePoints )
361 if ( !theFile.isOpen() )
364 // Strings in file is written as:
365 // 1. X(float) Y(float) Z(float)
366 // 2. X(float) Y(float) Z(float)
369 while ( !theFile.atEnd() )
371 QString aLine = theFile.readLine().simplified();
372 if ( aLine.isEmpty() )
375 QStringList aValues = aLine.split( QRegExp( "\\s+" ), QString::SkipEmptyParts );
376 if ( aValues.length() < 3 )
379 AltitudePoint aPoint;
381 QString anX = aValues.value( 0 );
382 QString anY = aValues.value( 1 );
383 QString aZ = aValues.value( 2 );
385 bool isXOk = false, isYOk = false, isZOk = false;
387 aPoint.SetX( anX.toDouble( &isXOk ) );
388 aPoint.SetY( anY.toDouble( &isYOk ) );
389 aPoint.SetZ( aZ.toDouble( &isZOk ) );
391 if ( !isXOk || !isYOk || !isZOk )