GRASS Programmer's Manual
6.4.2(2012)
|
00001 #include <rpc/types.h> 00002 #include <rpc/xdr.h> 00003 #include <stdio.h> 00004 #include <stdlib.h> 00005 #include <sys/types.h> 00006 #include <unistd.h> 00007 #include <grass/G3d.h> 00008 00009 /*--------------------------------------------------------------------------*/ 00010 00011 static unsigned char clearMask[9] = 00012 { 255, 128, 192, 224, 240, 248, 252, 254, 255 }; 00013 00014 /*---------------------------------------------------------------------------*/ 00015 00016 static void G3d_float2xdrFloat(float *f, float *xdrf) 00017 { 00018 XDR xdrEncodeStream; 00019 00020 xdrmem_create(&xdrEncodeStream, (caddr_t) xdrf, 4, XDR_ENCODE); 00021 00022 if (!xdr_setpos(&xdrEncodeStream, 0)) 00023 G3d_fatalError("G3d_float2xdrFloat: positioning xdr failed"); 00024 00025 if (!xdr_float(&xdrEncodeStream, f)) 00026 G3d_fatalError("G3d_float2xdrFloat: writing xdr failed"); 00027 00028 xdr_destroy(&xdrEncodeStream); 00029 } 00030 00031 /*---------------------------------------------------------------------------*/ 00032 00033 static void G3d_double2xdrDouble(double *d, double *xdrd) 00034 { 00035 XDR xdrEncodeStream; 00036 00037 xdrmem_create(&xdrEncodeStream, (caddr_t) xdrd, 8, XDR_ENCODE); 00038 00039 if (!xdr_setpos(&xdrEncodeStream, 0)) 00040 G3d_fatalError("G3d_double2xdrDouble: positioning xdr failed"); 00041 00042 if (!xdr_double(&xdrEncodeStream, d)) 00043 G3d_fatalError("G3d_double2xdrDouble: writing xdr failed"); 00044 00045 xdr_destroy(&xdrEncodeStream); 00046 } 00047 00048 /*---------------------------------------------------------------------------*/ 00049 00050 static void G3d_truncFloat(float *f, int p) 00051 { 00052 unsigned char *c; 00053 00054 if ((p == -1) || (p >= 23)) 00055 return; 00056 00057 c = (unsigned char *)f; 00058 00059 c++; 00060 if (p <= 7) { 00061 *c++ &= clearMask[(p + 1) % 8]; 00062 *c++ = 0; 00063 *c = 0; 00064 return; 00065 } 00066 00067 c++; 00068 if (p <= 15) { 00069 *c++ &= clearMask[(p + 1) % 8]; 00070 *c = 0; 00071 return; 00072 } 00073 00074 c++; 00075 *c &= clearMask[(p + 1) % 8]; 00076 return; 00077 } 00078 00079 /*---------------------------------------------------------------------------*/ 00080 00081 static void G3d_truncDouble(double *d, int p) 00082 { 00083 unsigned char *c; 00084 00085 if ((p == -1) || (p >= 52)) 00086 return; 00087 00088 c = (unsigned char *)d; 00089 00090 c++; 00091 if (p <= 4) { 00092 *c++ &= clearMask[(p + 4) % 8]; 00093 *c++ = 0; 00094 *c++ = 0; 00095 *c++ = 0; 00096 *c++ = 0; 00097 *c++ = 0; 00098 *c = 0; 00099 return; 00100 } 00101 00102 c++; 00103 if (p <= 12) { 00104 *c++ &= clearMask[(p + 4) % 8]; 00105 *c++ = 0; 00106 *c++ = 0; 00107 *c++ = 0; 00108 *c++ = 0; 00109 *c = 0; 00110 return; 00111 } 00112 00113 c++; 00114 if (p <= 20) { 00115 *c++ &= clearMask[(p + 4) % 8]; 00116 *c++ = 0; 00117 *c++ = 0; 00118 *c++ = 0; 00119 *c = 0; 00120 return; 00121 } 00122 00123 c++; 00124 if (p <= 28) { 00125 *c++ &= clearMask[(p + 4) % 8]; 00126 *c++ = 0; 00127 *c++ = 0; 00128 *c = 0; 00129 return; 00130 } 00131 00132 c++; 00133 if (p <= 36) { 00134 *c++ &= clearMask[(p + 4) % 8]; 00135 *c++ = 0; 00136 *c = 0; 00137 return; 00138 } 00139 00140 c++; 00141 if (p <= 44) { 00142 *c++ &= clearMask[(p + 4) % 8]; 00143 *c = 0; 00144 return; 00145 } 00146 00147 c++; 00148 *c &= clearMask[(p + 4) % 8]; 00149 return; 00150 } 00151 00152 /*---------------------------------------------------------------------------*/ 00153 00154 static void G3d_float2Double(float *f, double *d) 00155 { 00156 unsigned char *c1, *c2, sign, c; 00157 int e; 00158 00159 c1 = (unsigned char *)f; 00160 c2 = (unsigned char *)d; 00161 00162 sign = (*c1 & (unsigned char)128); 00163 e = (((*c1 & (unsigned char)127) << 1) | 00164 ((*(c1 + 1) & (unsigned char)128) >> 7)); 00165 00166 if ((*c1 != 0) || (*(c1 + 1) != 0) || (*(c1 + 2) != 0) || 00167 (*(c1 + 3) != 0)) 00168 e += 1023 - 127; 00169 c = e / 16; 00170 00171 *c2++ = (sign | c); 00172 00173 c1++; 00174 00175 c = e % 16; 00176 *c2 = (c << 4); 00177 *c2++ |= ((*c1 & (unsigned char)127) >> 3); 00178 00179 *c2 = ((*c1++ & (unsigned char)7) << 5); 00180 *c2++ |= (*c1 >> 3); 00181 00182 *c2 = ((*c1++ & (unsigned char)7) << 5); 00183 *c2++ |= (*c1 >> 3); 00184 00185 *c2++ = ((*c1 & (unsigned char)7) << 5); 00186 00187 *c2++ = (unsigned char)0; 00188 *c2++ = (unsigned char)0; 00189 *c2 = (unsigned char)0; 00190 } 00191 00192 /*---------------------------------------------------------------------------*/ 00193 00194 static int G3d_compareFloats(float *f1, int p1, float *f2, int p2) 00195 { 00196 unsigned char *c1, *c2; 00197 float xdrf1, xdrf2; 00198 00199 if (G3d_isNullValueNum(f1, FCELL_TYPE)) 00200 return G3d_isNullValueNum(f2, FCELL_TYPE); 00201 00202 G3d_float2xdrFloat(f1, &xdrf1); 00203 G3d_float2xdrFloat(f2, &xdrf2); 00204 00205 c1 = (unsigned char *)&xdrf1; 00206 c2 = (unsigned char *)&xdrf2; 00207 00208 /* printf ("%d %d (%d %d %d %d) (%d %d %d %d) %d\n", p1, p2, *c1, *(c1 + 1), *(c1 + 2), *(c1 + 3), *c2, *(c2 + 1), *(c2 + 2), *(c2 + 3), *f1 == *f2); */ 00209 00210 if ((p1 != -1) && (p1 < 23) && ((p1 < p2) || (p2 == -1))) 00211 G3d_truncFloat(&xdrf2, p1); 00212 if ((p2 != -1) && (p2 < 23) && ((p2 < p1) || (p1 == -1))) 00213 G3d_truncFloat(&xdrf1, p2); 00214 00215 /* printf ("%d %d (%d %d %d %d) (%d %d %d %d) %d\n", p1, p2, *c1, *(c1 + 1), *(c1 + 2), *(c1 + 3), *c2, *(c2 + 1), *(c2 + 2), *(c2 + 3), *f1 == *f2); */ 00216 00217 return (*c1 == *c2) && (*(c1 + 1) == *(c2 + 1)) && 00218 (*(c1 + 2) == *(c2 + 2)) 00219 && (*(c1 + 3) == *(c2 + 3)); 00220 } 00221 00222 00223 /*---------------------------------------------------------------------------*/ 00224 00225 static int G3d_compareDoubles(double *d1, int p1, double *d2, int p2) 00226 { 00227 unsigned char *c1, *c2; 00228 double xdrd1, xdrd2; 00229 00230 if (G3d_isNullValueNum(d1, DCELL_TYPE)) 00231 return G3d_isNullValueNum(d2, DCELL_TYPE); 00232 00233 G3d_double2xdrDouble(d1, &xdrd1); 00234 G3d_double2xdrDouble(d2, &xdrd2); 00235 00236 c1 = (unsigned char *)&xdrd1; 00237 c2 = (unsigned char *)&xdrd2; 00238 00239 /* printf ("%d %d (%d %d %d %d %d %d %d %d) (%d %d %d %d %d %d %d %d)\n", p1, p2, *c1, *(c1 + 1), *(c1 + 2), *(c1 + 3), *(c1 + 4), *(c1 + 5), *(c1 + 6), *(c1 + 7), *c2, *(c2 + 1), *(c2 + 2), *(c2 + 3), *(c2 + 4), *(c2 + 5), *(c2 + 6), *(c2 + 7)); */ 00240 00241 if ((p1 != -1) && (p1 < 52) && ((p1 < p2) || (p2 == -1))) 00242 G3d_truncDouble(&xdrd2, p1); 00243 if ((p2 != -1) && (p2 < 52) && ((p2 < p1) || (p1 == -1))) 00244 G3d_truncDouble(&xdrd1, p2); 00245 00246 /* printf ("%d %d (%d %d %d %d %d %d %d %d) (%d %d %d %d %d %d %d %d)\n", p1, p2, *c1, *(c1 + 1), *(c1 + 2), *(c1 + 3), *(c1 + 4), *(c1 + 5), *(c1 + 6), *(c1 + 7), *c2, *(c2 + 1), *(c2 + 2), *(c2 + 3), *(c2 + 4), *(c2 + 5), *(c2 + 6), *(c2 + 7)); */ 00247 00248 return (*c1 == *c2) && (*(c1 + 1) == *(c2 + 1)) && 00249 (*(c1 + 2) == *(c2 + 2)) 00250 && (*(c1 + 3) == *(c2 + 3)) && (*(c1 + 4) == *(c2 + 4)) 00251 && (*(c1 + 5) == *(c2 + 5)) && (*(c1 + 6) == *(c2 + 6)) 00252 && (*(c1 + 7) == *(c2 + 7)); 00253 } 00254 00255 00256 /*---------------------------------------------------------------------------*/ 00257 00258 static int G3d_compareFloatDouble(float *f, int p1, double *d, int p2) 00259 { 00260 unsigned char *c1, *c2; 00261 float xdrf, fTmp; 00262 double xdrd, xdrd2, dTmp; 00263 00264 if (G3d_isNullValueNum(f, FCELL_TYPE)) 00265 return G3d_isNullValueNum(d, DCELL_TYPE); 00266 00267 /* need this since assigning a double to a float actually may change the */ 00268 /* bit pattern. an example (in xdr format) is the double */ 00269 /* (63 237 133 81 81 108 3 32) which truncated to 23 bits precision should */ 00270 /* become (63 237 133 81 64 0 0 0). however assigned to a float it becomes */ 00271 /* (63 237 133 81 96 0 0 0). */ 00272 fTmp = *d; 00273 dTmp = fTmp; 00274 00275 G3d_float2xdrFloat(f, &xdrf); 00276 G3d_float2Double(&xdrf, &xdrd2); 00277 G3d_double2xdrDouble(&dTmp, &xdrd); 00278 00279 c1 = (unsigned char *)&xdrd2; 00280 c2 = (unsigned char *)&xdrd; 00281 00282 /* printf ("%d %d (%d %d %d %d) (%d %d %d %d %d %d %d %d) (%d %d %d %d %d %d %d %d)\n", p1, p2, *((unsigned char *) &xdrf), *(((unsigned char *) &xdrf) + 1), *(((unsigned char *) &xdrf) + 2), *(((unsigned char *) &xdrf) + 3), *c1, *(c1 + 1), *(c1 + 2), *(c1 + 3), *(c1 + 4), *(c1 + 5), *(c1 + 6), *(c1 + 7), *c2, *(c2 + 1), *(c2 + 2), *(c2 + 3), *(c2 + 4), *(c2 + 5), *(c2 + 6), *(c2 + 7)); */ 00283 00284 00285 if (((p1 != -1) && ((p1 < p2) || (p2 == -1))) || 00286 ((p1 == -1) && ((p2 > 23) || (p2 == -1)))) 00287 G3d_truncDouble(&xdrd, (p1 != -1 ? p1 : 23)); 00288 if ((p2 != -1) && (p2 < 23) && ((p2 < p1) || (p1 == -1))) 00289 G3d_truncDouble(&xdrd2, p2); 00290 00291 /* printf ("%d %d (%d %d %d %d) (%d %d %d %d %d %d %d %d) (%d %d %d %d %d %d %d %d)\n", p1, p2, *((unsigned char *) &xdrf), *(((unsigned char *) &xdrf) + 1), *(((unsigned char *) &xdrf) + 2), *(((unsigned char *) &xdrf) + 3), *c1, *(c1 + 1), *(c1 + 2), *(c1 + 3), *(c1 + 4), *(c1 + 5), *(c1 + 6), *(c1 + 7), *c2, *(c2 + 1), *(c2 + 2), *(c2 + 3), *(c2 + 4), *(c2 + 5), *(c2 + 6), *(c2 + 7)); */ 00292 00293 return (*c1 == *c2) && (*(c1 + 1) == *(c2 + 1)) && 00294 (*(c1 + 2) == *(c2 + 2)) 00295 && (*(c1 + 3) == *(c2 + 3)) && (*(c1 + 4) == *(c2 + 4)) 00296 && (*(c1 + 5) == *(c2 + 5)) && (*(c1 + 6) == *(c2 + 6)) 00297 && (*(c1 + 7) == *(c2 + 7)); 00298 } 00299 00300 /*---------------------------------------------------------------------------*/ 00301 00302 static void compareFilesNocache(void *map, void *map2) 00303 { 00304 double n1 = 0, n2 = 0; 00305 double *n1p, *n2p; 00306 float *f1p, *f2p; 00307 int x, y, z, correct; 00308 int p1, p2; 00309 int tileX, tileY, tileZ, typeIntern, typeIntern2; 00310 int nx, ny, nz; 00311 00312 p1 = G3d_tilePrecisionMap(map); 00313 p2 = G3d_tilePrecisionMap(map2); 00314 00315 G3d_getTileDimensionsMap(map, &tileX, &tileY, &tileZ); 00316 G3d_getNofTilesMap(map2, &nx, &ny, &nz); 00317 typeIntern = G3d_tileTypeMap(map); 00318 typeIntern2 = G3d_tileTypeMap(map2); 00319 00320 n1p = &n1; 00321 f1p = (float *)&n1; 00322 n2p = &n2; 00323 f2p = (float *)&n2; 00324 00325 for (z = 0; z < nz * tileZ; z++) { 00326 printf("comparing: z = %d\n", z); 00327 00328 for (y = 0; y < ny * tileY; y++) { 00329 for (x = 0; x < nx * tileX; x++) { 00330 00331 G3d_getBlock(map, x, y, z, 1, 1, 1, n1p, typeIntern); 00332 G3d_getBlock(map2, x, y, z, 1, 1, 1, n2p, typeIntern2); 00333 00334 if (typeIntern == FCELL_TYPE) { 00335 if (typeIntern2 == FCELL_TYPE) 00336 correct = G3d_compareFloats(f1p, p1, f2p, p2); 00337 else 00338 correct = G3d_compareFloatDouble(f1p, p1, n2p, p2); 00339 } 00340 else { 00341 if (typeIntern2 == FCELL_TYPE) 00342 correct = G3d_compareFloatDouble(f2p, p2, n1p, p1); 00343 else 00344 correct = G3d_compareDoubles(n1p, p1, n2p, p2); 00345 } 00346 00347 if (!correct) { 00348 int xTile, yTile, zTile, xOffs, yOffs, zOffs; 00349 00350 G3d_coord2tileCoord(map2, x, y, z, &xTile, &yTile, &zTile, 00351 &xOffs, &yOffs, &zOffs); 00352 printf("(%d %d %d) (%d %d %d) (%d %d %d) %.20f %.20f\n", 00353 x, y, z, xTile, yTile, zTile, xOffs, yOffs, zOffs, 00354 *n1p, *n2p); 00355 G3d_fatalError 00356 ("compareFilesNocache: files don't match\n"); 00357 } 00358 } 00359 } 00360 } 00361 00362 printf("Files are identical up to precision.\n"); 00363 } 00364 00365 /*---------------------------------------------------------------------------*/ 00366 00367 00385 void 00386 G3d_compareFiles(const char *f1, const char *mapset1, const char *f2, 00387 const char *mapset2) 00388 { 00389 void *map, *map2; 00390 double n1 = 0, n2 = 0; 00391 double *n1p, *n2p; 00392 float *f1p, *f2p; 00393 int x, y, z, correct; 00394 int p1, p2; 00395 int rows, cols, depths; 00396 int tileX, tileY, tileZ, typeIntern, typeIntern2, tileX2, tileY2, tileZ2; 00397 int nx, ny, nz; 00398 00399 printf("\nComparing %s and %s\n", f1, f2); 00400 00401 map = G3d_openCellOld(f1, mapset1, G3D_DEFAULT_WINDOW, 00402 G3D_TILE_SAME_AS_FILE, G3D_USE_CACHE_DEFAULT); 00403 if (map == NULL) 00404 G3d_fatalError("G3d_compareFiles: error in G3d_openCellOld"); 00405 00406 G3d_printHeader(map); 00407 00408 map2 = G3d_openCellOld(f2, mapset2, G3D_DEFAULT_WINDOW, 00409 G3D_TILE_SAME_AS_FILE, G3D_USE_CACHE_DEFAULT); 00410 if (map2 == NULL) 00411 G3d_fatalError("G3d_compareFiles: error in G3d_openCellOld"); 00412 00413 G3d_printHeader(map2); 00414 00415 typeIntern = G3d_tileTypeMap(map); 00416 typeIntern2 = G3d_tileTypeMap(map2); 00417 00418 p1 = G3d_tilePrecisionMap(map); 00419 p2 = G3d_tilePrecisionMap(map2); 00420 00421 G3d_getTileDimensionsMap(map, &tileX, &tileY, &tileZ); 00422 G3d_getTileDimensionsMap(map2, &tileX2, &tileY2, &tileZ2); 00423 G3d_getNofTilesMap(map2, &nx, &ny, &nz); 00424 G3d_getCoordsMap(map, &rows, &cols, &depths); 00425 00426 if ((!G3d_tileUseCacheMap(map)) || (!G3d_tileUseCacheMap(map2))) { 00427 compareFilesNocache(map, map2); 00428 G3d_closeCell(map); 00429 G3d_closeCell(map2); 00430 return; 00431 } 00432 00433 n1p = &n1; 00434 f1p = (float *)&n1; 00435 n2p = &n2; 00436 f2p = (float *)&n2; 00437 00438 G3d_autolockOn(map); 00439 G3d_autolockOn(map2); 00440 G3d_minUnlocked(map, cols / tileX + 1); 00441 00442 G3d_getCoordsMap(map2, &rows, &cols, &depths); 00443 G3d_minUnlocked(map2, cols / tileX + 1); 00444 00445 G3d_getCoordsMap(map, &rows, &cols, &depths); 00446 for (z = 0; z < depths; z++) { 00447 printf("comparing: z = %d\n", z); 00448 00449 if ((z % tileZ) == 0) { 00450 if (!G3d_unlockAll(map)) 00451 G3d_fatalError("G3d_compareFiles: error in G3d_unlockAll"); 00452 } 00453 if ((z % tileZ2) == 0) { 00454 if (!G3d_unlockAll(map2)) 00455 G3d_fatalError("G3d_compareFiles: error in G3d_unlockAll"); 00456 } 00457 00458 for (y = 0; y < rows; y++) { 00459 for (x = 0; x < cols; x++) { 00460 G3d_getValueRegion(map, x, y, z, n1p, typeIntern); 00461 G3d_getValueRegion(map2, x, y, z, n2p, typeIntern2); 00462 00463 G3d_isNullValueNum(n1p, typeIntern); 00464 G3d_isNullValueNum(n2p, typeIntern2); 00465 00466 if (typeIntern == FCELL_TYPE) { 00467 if (typeIntern2 == FCELL_TYPE) 00468 correct = G3d_compareFloats(f1p, p1, f2p, p2); 00469 else 00470 correct = G3d_compareFloatDouble(f1p, p1, n2p, p2); 00471 } 00472 else { 00473 if (typeIntern2 == FCELL_TYPE) 00474 correct = G3d_compareFloatDouble(f2p, p2, n1p, p1); 00475 else 00476 correct = G3d_compareDoubles(n1p, p1, n2p, p2); 00477 } 00478 00479 if (!correct) { 00480 int xTile, yTile, zTile, xOffs, yOffs, zOffs; 00481 00482 G3d_coord2tileCoord(map2, x, y, z, &xTile, &yTile, &zTile, 00483 &xOffs, &yOffs, &zOffs); 00484 G3d_fatalError("G3d_compareFiles: files don't match\n"); 00485 } 00486 } 00487 } 00488 } 00489 00490 printf("Files are identical up to precision.\n"); 00491 G3d_closeCell(map); 00492 G3d_closeCell(map2); 00493 }