[cig-commits] r12930 - cs/spatialdata-0.1/trunk/libsrc/spatialdb

[brad at geodynamics.org]

Author: brad
Date: 2008-09-21 12:57:33 -0700 (Sun, 21 Sep 2008)
New Revision: 12930

Modified:
   cs/spatialdata-0.1/trunk/libsrc/spatialdb/GocadVoxet.cc
   cs/spatialdata-0.1/trunk/libsrc/spatialdb/SCECCVMH.cc
Log:
Added retry loop to adjust search location in vertical direction when encountering 'holes' in the voxet file.

Modified: cs/spatialdata-0.1/trunk/libsrc/spatialdb/GocadVoxet.cc
===================================================================
--- cs/spatialdata-0.1/trunk/libsrc/spatialdb/GocadVoxet.cc	2008-09-21 11:33:29 UTC (rev 12929)
+++ cs/spatialdata-0.1/trunk/libsrc/spatialdb/GocadVoxet.cc	2008-09-21 19:57:33 UTC (rev 12930)
@@ -84,6 +84,22 @@
       indexZ >= 0 && indexZ < numZ) {
     int indexV = indexZ*numY*numX + indexY*numX + indexX;
     *value = _data[indexV];
+
+    // If voxet value is "no data value"
+    if (fabs(1.0 - *value / _property.noDataValue) < 1.0e-6) {
+      // If near indexZ=0, retry with indexZ+1, otherwise if near
+      // indexZ=numZ, retry with indexZ-1.
+      const int dz = (indexZ < numZ/2) ? +1 : -1;
+      const int maxRetries = 32;
+      for (int iTry=0; iTry < maxRetries; ++iTry) {
+	const int indexZNew = indexZ + dz*iTry;
+	assert(indexZNew >= 0 && indexZNew < numZ);
+	indexV = indexZNew*numY*numX + indexY*numX + indexX;
+	*value = _data[indexV];
+	if (fabs(1.0 - *value / _property.noDataValue) > 1.0e-6)
+	  break;
+      } // for
+    } // if
   } else {
     *value = _property.noDataValue;
     flag = 1;

Modified: cs/spatialdata-0.1/trunk/libsrc/spatialdb/SCECCVMH.cc
===================================================================
--- cs/spatialdata-0.1/trunk/libsrc/spatialdb/SCECCVMH.cc	2008-09-21 11:33:29 UTC (rev 12929)
+++ cs/spatialdata-0.1/trunk/libsrc/spatialdb/SCECCVMH.cc	2008-09-21 19:57:33 UTC (rev 12930)
@@ -235,21 +235,16 @@
       case QUERY_VP :
 	outsideVoxet = _queryVp(&vals[iVal]);
 	if (outsideVoxet)
-	  queryFlag = outsideVoxet || queryFlag;
+	  queryFlag |= outsideVoxet;
 	haveVp = true;
 	vp = vals[iVal];
 	break;
-      case QUERY_VPTAG :
-	outsideVoxet = _queryTag(&vals[iVal]);
-	if (outsideVoxet)
-	  queryFlag = outsideVoxet || queryFlag;
-	break;
       case QUERY_DENSITY :
 	if (!haveVp) {
 	  outsideVoxet = _queryVp(&vp);
 	  haveVp = true;
 	  if (outsideVoxet)
-	    queryFlag = outsideVoxet || queryFlag;
+	    queryFlag |= outsideVoxet;
 	} // if
 	vals[iVal] = _calcDensity(vp);
 	break;
@@ -258,7 +253,7 @@
 	  outsideVoxet = _queryVp(&vp);
 	  haveVp = true;
 	  if (outsideVoxet)
-	    queryFlag = outsideVoxet || queryFlag;
+	    queryFlag |= outsideVoxet;
 	} // if
 	vals[iVal] = _calcVs(vp);
 	break;
@@ -267,7 +262,7 @@
 	  assert(0 != _topoElev);
 	  outsideVoxet = _topoElev->queryNearest(&vals[iVal], _xyzUTM);
 	  if (outsideVoxet)
-	    queryFlag = outsideVoxet || queryFlag;
+	    queryFlag |= outsideVoxet;
 	} else
 	  vals[iVal] = topoElev;
 	break;
@@ -275,14 +270,19 @@
 	assert(0 != _baseDepth);
 	outsideVoxet = _baseDepth->queryNearest(&vals[iVal], _xyzUTM);
 	if (outsideVoxet)
-	  queryFlag = outsideVoxet || queryFlag;
+	  queryFlag |= outsideVoxet;
 	break;
       case QUERY_MOHODEPTH :
 	assert(0 != _mohoDepth);
 	outsideVoxet = _mohoDepth->queryNearest(&vals[iVal], _xyzUTM);
 	if (outsideVoxet)
-	  queryFlag = outsideVoxet || queryFlag;
+	  queryFlag |= outsideVoxet;
 	break;
+      case QUERY_VPTAG :
+	outsideVoxet = _queryTag(&vals[iVal]);
+	if (outsideVoxet)
+	  queryFlag |= outsideVoxet;
+	break;
       default:
 	assert(0);
       } // switch