Skip to content

Commit ee8ec43

Browse files
PauloCarvalhoRJPauloCarvalhoRJ
committed
ContactAnalysis: added support to GeoGrids.
1 parent fe9667e commit ee8ec43

File tree

1 file changed

+55
-23
lines changed

1 file changed

+55
-23
lines changed

‎geostats/contactanalysis.cpp

Lines changed: 55 additions & 23 deletions
Original file line numberDiff line numberDiff line change
@@ -9,12 +9,13 @@
99
#include "domain/segmentset.h"
1010
#include "domain/cartesiangrid.h"
1111
#include "domain/segmentset.h"
12+
#include "domain/geogrid.h"
1213
#include "spatialindex/spatialindex.h"
1314
#include "geostats/searchannulus.h"
1415
#include "geostats/searchwasher.h"
1516
#include "geostats/searchverticaldumbbell.h"
1617
#include "geostats/searchstrategy.h"
17-
#include "geostats/searchhollowsphere.h"
18+
#include "geostats/searchsphericalshell.h"
1819
#include "geostats/pointsetcell.h"
1920
#include "geostats/gridcell.h"
2021
#include "geostats/segmentsetcell.h"
@@ -175,35 +176,40 @@ DataCellPtrMultiset ContactAnalysis::getSamplesFromInputDataSet(const DataCell &
175176
for( ; it != samplesIndexes.end(); ++it ){
176177
switch ( m_inputDataType ) {
177178
case InputDataSetType::POINTSET:
178-
{
179+
{
179180
DataCellPtr p(new PointSetCell( static_cast<PointSet*>( m_inputDataFile ), m_attributeGrade->getAttributeGEOEASgivenIndex()-1, *it ));
180181
p->computeCartesianDistance( sample );
181182
result.insert( p );
182-
}
183+
}
183184
break;
184185
case InputDataSetType::CARTESIANGRID:
185-
{
186+
{
186187
CartesianGrid* cg = static_cast<CartesianGrid*>( m_inputDataFile );
187188
uint i, j, k;
188189
cg->indexToIJK( *it, i, j, k );
189190
DataCellPtr p(new GridCell( cg, m_attributeGrade->getAttributeGEOEASgivenIndex()-1, i, j, k ));
190191
p->computeCartesianDistance( sample );
191192
result.insert( p );
192-
}
193+
}
193194
break;
194195
case InputDataSetType::GEOGRID:
195-
{
196-
Application::instance()->logError( "ContactAnalysis::getSamplesFromInputDataSet(): GeoGrids cannot be used as input data yet."
197-
" It is necessary to create a class like GeoGridCell inheriting from DataCell." );
198-
return result;
199-
}
196+
{
197+
GeoGrid* ggAspect = dynamic_cast<GeoGrid*>(m_inputDataFile);
198+
uint i, j, k;
199+
ggAspect->indexToIJK( *it, i, j, k );
200+
DataCellPtr p( new GridCell( ggAspect,
201+
m_attributeGrade->getAttributeGEOEASgivenIndex()-1,
202+
i, j, k ) );
203+
p->computeCartesianDistance( sample );
204+
result.insert( p );
205+
}
200206
break;
201207
case InputDataSetType::SEGMENTSET:
202-
{
208+
{
203209
DataCellPtr p(new SegmentSetCell( static_cast<SegmentSet*>( m_inputDataFile ), m_attributeGrade->getAttributeGEOEASgivenIndex()-1, *it ));
204210
p->computeCartesianDistance( sample );
205211
result.insert( p );
206-
}
212+
}
207213
break;
208214
default:
209215
Application::instance()->logError( "ContactAnalysis::getSamplesFromInputDataSet(): Input data file type not recognized or undefined." );
@@ -284,6 +290,9 @@ bool ContactAnalysis::run()
284290
} else if ( m_inputDataFile->getFileType() == "CARTESIANGRID") {
285291
CartesianGrid* cgAspect = dynamic_cast<CartesianGrid*>( m_inputDataFile );
286292
spatialIndex.fill( cgAspect );
293+
} else if ( m_inputDataFile->getFileType() == "GEOGRID") {
294+
GeoGrid* ggAspect = dynamic_cast<GeoGrid*>( m_inputDataFile );
295+
spatialIndex.fillWithCenters( ggAspect, 0.0001 );
287296
} else {
288297
m_lastError = "Internal error building spatial index: input data of type " + m_inputDataFile->getFileType() + " are not currently supported.";
289298
return false;
@@ -349,8 +358,30 @@ bool ContactAnalysis::run()
349358
m_lastError = "Cannot perform vertical contact analysis on a 2D dataset.";
350359
return false;
351360
} else { //if data set is 3D
352-
if( ! m_inputDataFile->isGridded() || m_inputDataType == InputDataSetType::CARTESIANGRID ) {
353-
searchNeighborhood.reset( new SearchVerticalDumbbell( m_lagSize, 2*current_lag, m_lagSize ) );
361+
if( ! m_inputDataFile->isGridded() ) {
362+
searchNeighborhood.reset( new SearchVerticalDumbbell( m_lagSize, 2*(current_lag-m_lagSize), m_lagSize ) );
363+
} else if ( m_inputDataType == InputDataSetType::CARTESIANGRID ) {
364+
CartesianGrid* cg = dynamic_cast< CartesianGrid* >( m_inputDataFile );
365+
const double dx = cg->getDX();
366+
const double dy = cg->getDY();
367+
double radius = std::sqrt( dx*dx + dy*dy ) / 2;
368+
searchNeighborhood.reset( new SearchVerticalDumbbell( m_lagSize, 2*(current_lag-m_lagSize), radius ) );
369+
} else if ( m_inputDataType == InputDataSetType::GEOGRID ) {
370+
GeoGrid* gg = dynamic_cast< GeoGrid* >( m_inputDataFile );
371+
double bbMinX, bbMinY, bbMinZ; //BB == Bounding Box
372+
double bbMaxX, bbMaxY, bbMaxZ;
373+
double maxBBdX = -std::numeric_limits<double>::max();
374+
double maxBBdY = -std::numeric_limits<double>::max();
375+
uint maxIndex = gg->getDataLineCount();
376+
for( uint i = 0; i < maxIndex; i++ ){
377+
gg->getBoundingBox( i, bbMinX, bbMinY, bbMinZ, bbMaxX, bbMaxY, bbMaxZ );
378+
double bbDx = bbMaxX - bbMinX;
379+
double bbDy = bbMaxY - bbMinY;
380+
if( bbDx > maxBBdX ) maxBBdX = bbDx;
381+
if( bbDy > maxBBdY ) maxBBdY = bbDy;
382+
}
383+
double radius = std::sqrt( maxBBdX*maxBBdX + maxBBdY*maxBBdY ) / 2;
384+
searchNeighborhood.reset( new SearchVerticalDumbbell( m_lagSize, 2*(current_lag-m_lagSize), radius ) );
354385
} else {
355386
m_lastError = "Internal error: no search neighborhood is available "
356387
"for vertical contact analysis with datasets of type " + m_inputDataFile->getFileType() + ".";
@@ -362,13 +393,7 @@ bool ContactAnalysis::run()
362393
m_lastError = "Cannot perform omnidirectional contact analysis on a 2D dataset.";
363394
return false;
364395
} else { //if data set is 3D
365-
if( ! m_inputDataFile->isGridded() || m_inputDataType == InputDataSetType::CARTESIANGRID ) {
366-
searchNeighborhood.reset( new SearchHollowSphere( current_lag - m_lagSize, current_lag ) );
367-
} else {
368-
m_lastError = "Internal error: no search neighborhood is available "
369-
"for omnidirectional contact analysis with datasets of type " + m_inputDataFile->getFileType() + ".";
370-
return false;
371-
}
396+
searchNeighborhood.reset( new SearchSphericalShell( current_lag - m_lagSize, current_lag ) );
372397
}
373398
}
374399
}
@@ -412,8 +437,15 @@ bool ContactAnalysis::run()
412437
}
413438
break;
414439
case InputDataSetType::GEOGRID:
415-
m_lastError = "Internal error: geogrids are not currently supported in contact analysis.";
416-
return false;
440+
{
441+
GeoGrid* ggAspect = dynamic_cast<GeoGrid*>(m_inputDataFile);
442+
uint i, j, k;
443+
ggAspect->indexToIJK( iCurrentSample, i, j, k );
444+
currentSampleCell.reset( new GridCell( ggAspect,
445+
m_attributeGrade->getAttributeGEOEASgivenIndex()-1,
446+
i, j, k ) );
447+
}
448+
break;
417449
case InputDataSetType::SEGMENTSET:
418450
currentSampleCell.reset( new SegmentSetCell( dynamic_cast<SegmentSet*>(m_inputDataFile),
419451
m_attributeGrade->getAttributeGEOEASgivenIndex()-1,

0 commit comments

Comments
 (0)