00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018 #include "qgsinterpolator.h"
00019 #include "qgsvectordataprovider.h"
00020 #include "qgsgeometry.h"
00021 #ifndef Q_OS_MACX
00022 #include <cmath>
00023 #else
00024 #include <math.h>
00025 #endif
00026 #ifdef WIN32
00027 #include <float.h>
00028 #define isnan(f) _isnan(f)
00029 #endif
00030
00031 QgsInterpolator::QgsInterpolator( const QList<LayerData>& layerData ): mDataIsCached( false ), mLayerData( layerData )
00032 {
00033
00034 }
00035
00036 QgsInterpolator::QgsInterpolator()
00037 {
00038
00039 }
00040
00041 QgsInterpolator::~QgsInterpolator()
00042 {
00043
00044 }
00045
00046 int QgsInterpolator::cacheBaseData()
00047 {
00048 if ( mLayerData.size() < 1 )
00049 {
00050 return 0;
00051 }
00052
00053
00054 mCachedBaseData.clear();
00055 mCachedBaseData.reserve( 100000 );
00056
00057 QList<LayerData>::iterator v_it = mLayerData.begin();
00058
00059 for ( ; v_it != mLayerData.end(); ++v_it )
00060 {
00061 if ( v_it->vectorLayer == 0 )
00062 {
00063 continue;
00064 }
00065
00066 QgsVectorDataProvider* provider = v_it->vectorLayer->dataProvider();
00067 if ( !provider )
00068 {
00069 return 2;
00070 }
00071
00072 QgsAttributeList attList;
00073 if ( !v_it->zCoordInterpolation )
00074 {
00075 attList.push_back( v_it->interpolationAttribute );
00076 }
00077
00078 provider->select( attList );
00079
00080 QgsFeature theFeature;
00081 double attributeValue = 0.0;
00082 bool attributeConversionOk = false;
00083
00084 while ( provider->nextFeature( theFeature ) )
00085 {
00086 if ( !v_it->zCoordInterpolation )
00087 {
00088 QgsAttributeMap attMap = theFeature.attributeMap();
00089 QgsAttributeMap::const_iterator att_it = attMap.find( v_it->interpolationAttribute );
00090 if ( att_it == attMap.end() )
00091 {
00092 continue;
00093 }
00094 attributeValue = att_it.value().toDouble( &attributeConversionOk );
00095 if ( !attributeConversionOk || isnan( attributeValue ) )
00096 {
00097 continue;
00098 }
00099 }
00100
00101 if ( addVerticesToCache( theFeature.geometry(), v_it->zCoordInterpolation, attributeValue ) != 0 )
00102 {
00103 return 3;
00104 }
00105 }
00106 }
00107
00108 return 0;
00109 }
00110
00111 int QgsInterpolator::addVerticesToCache( QgsGeometry* geom, bool zCoord, double attributeValue )
00112 {
00113 if ( !geom )
00114 {
00115 return 1;
00116 }
00117
00118 bool hasZValue = false;
00119 unsigned char* currentWkbPtr = geom->asWkb();
00120 vertexData theVertex;
00121
00122 QGis::WkbType wkbType = geom->wkbType();
00123 switch ( wkbType )
00124 {
00125 case QGis::WKBPoint25D:
00126 hasZValue = true;
00127 case QGis::WKBPoint:
00128 {
00129 currentWkbPtr += ( 1 + sizeof( int ) );
00130 theVertex.x = *(( double * )( currentWkbPtr ) );
00131 currentWkbPtr += sizeof( double );
00132 theVertex.y = *(( double * )( currentWkbPtr ) );
00133 if ( zCoord && hasZValue )
00134 {
00135 currentWkbPtr += sizeof( double );
00136 theVertex.z = *(( double * )( currentWkbPtr ) );
00137 }
00138 else
00139 {
00140 theVertex.z = attributeValue;
00141 }
00142 mCachedBaseData.push_back( theVertex );
00143 break;
00144 }
00145 case QGis::WKBLineString25D:
00146 hasZValue = true;
00147 case QGis::WKBLineString:
00148 {
00149 currentWkbPtr += ( 1 + sizeof( int ) );
00150 int* npoints = ( int* )currentWkbPtr;
00151 currentWkbPtr += sizeof( int );
00152 for ( int index = 0; index < *npoints; ++index )
00153 {
00154 theVertex.x = *(( double * )( currentWkbPtr ) );
00155 currentWkbPtr += sizeof( double );
00156 theVertex.y = *(( double * )( currentWkbPtr ) );
00157 currentWkbPtr += sizeof( double );
00158 if ( zCoord && hasZValue )
00159 {
00160 theVertex.z = *(( double * )( currentWkbPtr ) );
00161 currentWkbPtr += sizeof( double );
00162 }
00163 else
00164 {
00165 theVertex.z = attributeValue;
00166 }
00167 mCachedBaseData.push_back( theVertex );
00168 }
00169 break;
00170 }
00171 #if 0
00172 case QGis::WKBPolygon25D:
00173 hasZValue = true;
00174 case QGis::WKBPolygon:
00175 {
00176 int* nrings = ( int* )( mGeometry + 5 );
00177 int* npoints;
00178 unsigned char* ptr = mGeometry + 9;
00179 for ( int index = 0; index < *nrings; ++index )
00180 {
00181 npoints = ( int* )ptr;
00182 ptr += sizeof( int );
00183 for ( int index2 = 0; index2 < *npoints; ++index2 )
00184 {
00185 tempx = ( double* )ptr;
00186 ptr += sizeof( double );
00187 tempy = ( double* )ptr;
00188 if ( point.sqrDist( *tempx, *tempy ) < actdist )
00189 {
00190 x = *tempx;
00191 y = *tempy;
00192 actdist = point.sqrDist( *tempx, *tempy );
00193 vertexnr = vertexcounter;
00194
00195 if ( index2 == 0 )
00196 {
00197 beforeVertex = vertexcounter + ( *npoints - 2 );
00198 afterVertex = vertexcounter + 1;
00199 }
00200 else if ( index2 == ( *npoints - 1 ) )
00201 {
00202 beforeVertex = vertexcounter - 1;
00203 afterVertex = vertexcounter - ( *npoints - 2 );
00204 }
00205 else
00206 {
00207 beforeVertex = vertexcounter - 1;
00208 afterVertex = vertexcounter + 1;
00209 }
00210 }
00211 ptr += sizeof( double );
00212 if ( hasZValue )
00213 {
00214 ptr += sizeof( double );
00215 }
00216 ++vertexcounter;
00217 }
00218 }
00219 break;
00220 }
00221 case QGis::WKBMultiPoint25D:
00222 hasZValue = true;
00223 case QGis::WKBMultiPoint:
00224 {
00225 unsigned char* ptr = mGeometry + 5;
00226 int* npoints = ( int* )ptr;
00227 ptr += sizeof( int );
00228 for ( int index = 0; index < *npoints; ++index )
00229 {
00230 ptr += ( 1 + sizeof( int ) );
00231 tempx = ( double* )ptr;
00232 tempy = ( double* )( ptr + sizeof( double ) );
00233 if ( point.sqrDist( *tempx, *tempy ) < actdist )
00234 {
00235 x = *tempx;
00236 y = *tempy;
00237 actdist = point.sqrDist( *tempx, *tempy );
00238 vertexnr = index;
00239 }
00240 ptr += ( 2 * sizeof( double ) );
00241 if ( hasZValue )
00242 {
00243 ptr += sizeof( double );
00244 }
00245 }
00246 break;
00247 }
00248 case QGis::WKBMultiLineString25D:
00249 hasZValue = true;
00250 case QGis::WKBMultiLineString:
00251 {
00252 unsigned char* ptr = mGeometry + 5;
00253 int* nlines = ( int* )ptr;
00254 int* npoints = 0;
00255 ptr += sizeof( int );
00256 for ( int index = 0; index < *nlines; ++index )
00257 {
00258 ptr += ( sizeof( int ) + 1 );
00259 npoints = ( int* )ptr;
00260 ptr += sizeof( int );
00261 for ( int index2 = 0; index2 < *npoints; ++index2 )
00262 {
00263 tempx = ( double* )ptr;
00264 ptr += sizeof( double );
00265 tempy = ( double* )ptr;
00266 ptr += sizeof( double );
00267 if ( point.sqrDist( *tempx, *tempy ) < actdist )
00268 {
00269 x = *tempx;
00270 y = *tempy;
00271 actdist = point.sqrDist( *tempx, *tempy );
00272 vertexnr = vertexcounter;
00273
00274 if ( index2 == 0 )
00275 {
00276 beforeVertex = -1;
00277 }
00278 else
00279 {
00280 beforeVertex = vertexnr - 1;
00281 }
00282 if ( index2 == ( *npoints ) - 1 )
00283 {
00284 afterVertex = -1;
00285 }
00286 else
00287 {
00288 afterVertex = vertexnr + 1;
00289 }
00290 }
00291 if ( hasZValue )
00292 {
00293 ptr += sizeof( double );
00294 }
00295 ++vertexcounter;
00296 }
00297 }
00298 break;
00299 }
00300 case QGis::WKBMultiPolygon25D:
00301 hasZValue = true;
00302 case QGis::WKBMultiPolygon:
00303 {
00304 unsigned char* ptr = mGeometry + 5;
00305 int* npolys = ( int* )ptr;
00306 int* nrings;
00307 int* npoints;
00308 ptr += sizeof( int );
00309 for ( int index = 0; index < *npolys; ++index )
00310 {
00311 ptr += ( 1 + sizeof( int ) );
00312 nrings = ( int* )ptr;
00313 ptr += sizeof( int );
00314 for ( int index2 = 0; index2 < *nrings; ++index2 )
00315 {
00316 npoints = ( int* )ptr;
00317 ptr += sizeof( int );
00318 for ( int index3 = 0; index3 < *npoints; ++index3 )
00319 {
00320 tempx = ( double* )ptr;
00321 ptr += sizeof( double );
00322 tempy = ( double* )ptr;
00323 if ( point.sqrDist( *tempx, *tempy ) < actdist )
00324 {
00325 x = *tempx;
00326 y = *tempy;
00327 actdist = point.sqrDist( *tempx, *tempy );
00328 vertexnr = vertexcounter;
00329
00330
00331 if ( index3 == 0 )
00332 {
00333 beforeVertex = vertexcounter + ( *npoints - 2 );
00334 afterVertex = vertexcounter + 1;
00335 }
00336 else if ( index3 == ( *npoints - 1 ) )
00337 {
00338 beforeVertex = vertexcounter - 1;
00339 afterVertex = vertexcounter - ( *npoints - 2 );
00340 }
00341 else
00342 {
00343 beforeVertex = vertexcounter - 1;
00344 afterVertex = vertexcounter + 1;
00345 }
00346 }
00347 ptr += sizeof( double );
00348 if ( hasZValue )
00349 {
00350 ptr += sizeof( double );
00351 }
00352 ++vertexcounter;
00353 }
00354 }
00355 }
00356 break;
00357 }
00358 #endif //0
00359 default:
00360 break;
00361 }
00362 mDataIsCached = true;
00363 return 0;
00364 }