GRASS Programmer's Manual 6.4.1(2011)
|
00001 00016 #include <string.h> 00017 #include <unistd.h> 00018 #include <sys/types.h> 00019 00020 #include <rpc/types.h> /* need this for sgi */ 00021 #include <rpc/xdr.h> 00022 00023 #include <grass/config.h> 00024 #include <grass/glocale.h> 00025 00026 #include "G.h" 00027 00028 /*--------------------------------------------------------------------------*/ 00029 00030 #define NULL_FILE "null" 00031 00032 /*--------------------------------------------------------------------------*/ 00033 00034 static int embed_nulls(int, void *, int, RASTER_MAP_TYPE, int, int); 00035 00036 /*--------------------------------------------------------------------------*/ 00037 00038 static int compute_window_row(int fd, int row, int *cellRow) 00039 { 00040 struct fileinfo *fcb = &G__.fileinfo[fd]; 00041 double f; 00042 int r; 00043 00044 /* check for row in window */ 00045 if (row < 0 || row >= G__.window.rows) { 00046 G_warning(_("Reading raster map <%s@%s> request for row %d is outside region"), 00047 fcb->name, fcb->mapset, row); 00048 00049 return -1; 00050 } 00051 00052 /* convert window row to cell file row */ 00053 f = row * fcb->C1 + fcb->C2; 00054 r = (int)f; 00055 if (f < r) /* adjust for rounding up of negatives */ 00056 r--; 00057 00058 if (r < 0 || r >= fcb->cellhd.rows) 00059 return 0; 00060 00061 *cellRow = r; 00062 00063 return 1; 00064 } 00065 00066 /*--------------------------------------------------------------------------*/ 00067 00068 static void do_reclass_int(int fd, void *cell, int null_is_zero) 00069 { 00070 struct fileinfo *fcb = &G__.fileinfo[fd]; 00071 CELL *c = cell; 00072 CELL *reclass_table = fcb->reclass.table; 00073 CELL min = fcb->reclass.min; 00074 CELL max = fcb->reclass.max; 00075 int i; 00076 00077 for (i = 0; i < G__.window.cols; i++) { 00078 if (G_is_c_null_value(&c[i])) { 00079 if (null_is_zero) 00080 c[i] = 0; 00081 continue; 00082 } 00083 00084 if (c[i] < min || c[i] > max) { 00085 if (null_is_zero) 00086 c[i] = 0; 00087 else 00088 G_set_c_null_value(&c[i], 1); 00089 continue; 00090 } 00091 00092 c[i] = reclass_table[c[i] - min]; 00093 00094 if (null_is_zero && G_is_c_null_value(&c[i])) 00095 c[i] = 0; 00096 } 00097 } 00098 00099 /*--------------------------------------------------------------------------*/ 00100 00101 static int read_data_fp_compressed(int fd, int row, unsigned char *data_buf, 00102 int *nbytes) 00103 { 00104 struct fileinfo *fcb = &G__.fileinfo[fd]; 00105 off_t t1 = fcb->row_ptr[row]; 00106 off_t t2 = fcb->row_ptr[row + 1]; 00107 size_t readamount = t2 - t1; 00108 size_t bufsize = fcb->cellhd.cols * fcb->nbytes; 00109 00110 if (lseek(fd, t1, SEEK_SET) < 0) 00111 return -1; 00112 00113 *nbytes = fcb->nbytes; 00114 00115 if ((size_t) G_zlib_read(fd, readamount, data_buf, bufsize) != bufsize) 00116 return -1; 00117 00118 return 0; 00119 } 00120 00121 /*--------------------------------------------------------------------------*/ 00122 00123 static void rle_decompress(unsigned char *dst, const unsigned char *src, 00124 int nbytes, int size) 00125 { 00126 int pairs = size / (nbytes + 1); 00127 int i; 00128 00129 for (i = 0; i < pairs; i++) { 00130 int repeat = *src++; 00131 int j; 00132 00133 for (j = 0; j < repeat; j++) { 00134 memcpy(dst, src, nbytes); 00135 dst += nbytes; 00136 } 00137 00138 src += nbytes; 00139 } 00140 } 00141 00142 static int read_data_compressed(int fd, int row, unsigned char *data_buf, 00143 int *nbytes) 00144 { 00145 struct fileinfo *fcb = &G__.fileinfo[fd]; 00146 off_t t1 = fcb->row_ptr[row]; 00147 off_t t2 = fcb->row_ptr[row + 1]; 00148 ssize_t readamount = t2 - t1; 00149 unsigned char *cmp = G__.compressed_buf; 00150 int n; 00151 00152 if (lseek(fd, t1, SEEK_SET) < 0) 00153 return -1; 00154 00155 if (read(fd, cmp, readamount) != readamount) 00156 return -1; 00157 00158 /* Now decompress the row */ 00159 if (fcb->cellhd.compressed > 0) { 00160 /* one byte is nbyte count */ 00161 n = *nbytes = *cmp++; 00162 readamount--; 00163 } 00164 else 00165 /* pre 3.0 compression */ 00166 n = *nbytes = fcb->nbytes; 00167 00168 if (fcb->cellhd.compressed < 0 || readamount < n * fcb->cellhd.cols) { 00169 if (fcb->cellhd.compressed == 2) 00170 G_zlib_expand(cmp, readamount, data_buf, n * fcb->cellhd.cols); 00171 else 00172 rle_decompress(data_buf, cmp, n, readamount); 00173 } 00174 else 00175 memcpy(data_buf, cmp, readamount); 00176 00177 return 0; 00178 } 00179 00180 /*--------------------------------------------------------------------------*/ 00181 00182 static int read_data_uncompressed(int fd, int row, unsigned char *data_buf, 00183 int *nbytes) 00184 { 00185 struct fileinfo *fcb = &G__.fileinfo[fd]; 00186 ssize_t bufsize = fcb->cellhd.cols * fcb->nbytes; 00187 00188 *nbytes = fcb->nbytes; 00189 00190 if (lseek(fd, (off_t) row * bufsize, SEEK_SET) == -1) 00191 return -1; 00192 00193 if (read(fd, data_buf, bufsize) != bufsize) 00194 return -1; 00195 00196 return 0; 00197 } 00198 00199 /*--------------------------------------------------------------------------*/ 00200 00201 #ifdef HAVE_GDAL 00202 static int read_data_gdal(int fd, int row, unsigned char *data_buf, int *nbytes) 00203 { 00204 struct fileinfo *fcb = &G__.fileinfo[fd]; 00205 CPLErr err; 00206 00207 *nbytes = fcb->nbytes; 00208 00209 err = G_gdal_raster_IO( 00210 fcb->gdal->band, GF_Read, 0, row, fcb->cellhd.cols, 1, data_buf, 00211 fcb->cellhd.cols, 1, fcb->gdal->type, 0, 0); 00212 00213 return err == CE_None ? 0 : -1; 00214 } 00215 #endif 00216 00217 /*--------------------------------------------------------------------------*/ 00218 00219 /* Actually read a row of data in */ 00220 00221 static int read_data(int fd, int row, unsigned char *data_buf, int *nbytes) 00222 { 00223 struct fileinfo *fcb = &G__.fileinfo[fd]; 00224 00225 #ifdef HAVE_GDAL 00226 if (fcb->gdal) 00227 return read_data_gdal(fd, row, data_buf, nbytes); 00228 #endif 00229 00230 if (!fcb->cellhd.compressed) 00231 return read_data_uncompressed(fd, row, data_buf, nbytes); 00232 00233 /* map is in compressed form */ 00234 00235 if (fcb->map_type == CELL_TYPE) 00236 return read_data_compressed(fd, row, data_buf, nbytes); 00237 else 00238 return read_data_fp_compressed(fd, row, data_buf, nbytes); 00239 } 00240 00241 /*--------------------------------------------------------------------------*/ 00242 00243 /* copy cell file data to user buffer translated by window column mapping */ 00244 00245 static void cell_values_int(int fd, const unsigned char *data, 00246 const COLUMN_MAPPING * cmap, int nbytes, 00247 void *cell, int n) 00248 { 00249 CELL *c = cell; 00250 COLUMN_MAPPING cmapold = 0; 00251 int big = (size_t) nbytes >= sizeof(CELL); 00252 int i; 00253 00254 for (i = 0; i < n; i++) { 00255 const unsigned char *d; 00256 int neg; 00257 CELL v; 00258 int j; 00259 00260 if (!cmap[i]) { 00261 c[i] = 0; 00262 continue; 00263 } 00264 00265 if (cmap[i] == cmapold) { 00266 c[i] = c[i - 1]; 00267 continue; 00268 } 00269 00270 d = data + (cmap[i] - 1) * nbytes; 00271 00272 if (big && (*d & 0x80)) { 00273 neg = 1; 00274 v = *d++ & 0x7f; 00275 } 00276 else { 00277 neg = 0; 00278 v = *d++; 00279 } 00280 00281 for (j = 1; j < nbytes; j++) 00282 v = (v << 8) + *d++; 00283 00284 c[i] = neg ? -v : v; 00285 00286 cmapold = cmap[i]; 00287 } 00288 } 00289 00290 /*--------------------------------------------------------------------------*/ 00291 00292 static void cell_values_float(int fd, const unsigned char *data, 00293 const COLUMN_MAPPING * cmap, int nbytes, 00294 void *cell, int n) 00295 { 00296 struct fileinfo *fcb = &G__.fileinfo[fd]; 00297 FCELL *c = cell; 00298 COLUMN_MAPPING cmapold = 0; 00299 XDR *xdrs = &fcb->xdrstream; 00300 int i; 00301 00302 /* xdr stream is initialized to read from */ 00303 /* fcb->data in 'opencell.c' */ 00304 xdr_setpos(xdrs, 0); 00305 00306 for (i = 0; i < n; i++) { 00307 if (!cmap[i]) { 00308 c[i] = 0; 00309 continue; 00310 } 00311 00312 if (cmap[i] == cmapold) { 00313 c[i] = c[i - 1]; 00314 continue; 00315 } 00316 00317 if (cmap[i] < cmapold) { 00318 xdr_setpos(xdrs, 0); 00319 cmapold = 0; 00320 } 00321 00322 while (cmapold++ != cmap[i]) /* skip */ 00323 if (!xdr_float(xdrs, &c[i])) 00324 G_fatal_error(_("cell_values_float: xdr_float failed for index %d"), 00325 i); 00326 00327 cmapold--; 00328 } 00329 } 00330 00331 /*--------------------------------------------------------------------------*/ 00332 00333 static void cell_values_double(int fd, const unsigned char *data, 00334 const COLUMN_MAPPING * cmap, int nbytes, 00335 void *cell, int n) 00336 { 00337 struct fileinfo *fcb = &G__.fileinfo[fd]; 00338 DCELL *c = cell; 00339 COLUMN_MAPPING cmapold = 0; 00340 XDR *xdrs = &fcb->xdrstream; 00341 int i; 00342 00343 /* xdr stream is initialized to read from */ 00344 /* fcb->data in 'opencell.c' */ 00345 xdr_setpos(xdrs, 0); 00346 00347 for (i = 0; i < n; i++) { 00348 if (!cmap[i]) { 00349 c[i] = 0; 00350 continue; 00351 } 00352 00353 if (cmap[i] == cmapold) { 00354 c[i] = c[i - 1]; 00355 continue; 00356 } 00357 00358 if (cmap[i] < cmapold) { 00359 xdr_setpos(xdrs, 0); 00360 cmapold = 0; 00361 } 00362 00363 while (cmapold++ != cmap[i]) /* skip */ 00364 if (!xdr_double(xdrs, &c[i])) 00365 G_fatal_error(_("cell_values_double: xdr_double failed for index %d"), 00366 i); 00367 00368 cmapold--; 00369 } 00370 } 00371 00372 /*--------------------------------------------------------------------------*/ 00373 00374 /*--------------------------------------------------------------------------*/ 00375 00376 #ifdef HAVE_GDAL 00377 00378 /*--------------------------------------------------------------------------*/ 00379 00380 static void gdal_values_int(int fd, const unsigned char *data, 00381 const COLUMN_MAPPING *cmap, int nbytes, 00382 CELL *cell, int n) 00383 { 00384 struct fileinfo *fcb = &G__.fileinfo[fd]; 00385 const unsigned char *d; 00386 COLUMN_MAPPING cmapold = 0; 00387 int i; 00388 00389 for (i = 0; i < n; i++) { 00390 if (!cmap[i]) { 00391 cell[i] = 0; 00392 continue; 00393 } 00394 00395 if (cmap[i] == cmapold) { 00396 cell[i] = cell[i-1]; 00397 continue; 00398 } 00399 00400 d = data + (cmap[i] - 1) * nbytes; 00401 00402 switch (fcb->gdal->type) { 00403 case GDT_Byte: cell[i] = *(GByte *)d; break; 00404 case GDT_Int16: cell[i] = *(GInt16 *)d; break; 00405 case GDT_UInt16: cell[i] = *(GUInt16 *)d; break; 00406 case GDT_Int32: cell[i] = *(GInt32 *)d; break; 00407 case GDT_UInt32: cell[i] = *(GUInt32 *)d; break; 00408 default: 00409 /* shouldn't happen */ 00410 G_set_c_null_value(&cell[i], 1); 00411 break; 00412 } 00413 00414 cmapold = cmap[i]; 00415 } 00416 } 00417 00418 /*--------------------------------------------------------------------------*/ 00419 00420 static void gdal_values_float(int fd, const float *data, 00421 const COLUMN_MAPPING *cmap, int nbytes, 00422 FCELL *cell, int n) 00423 { 00424 COLUMN_MAPPING cmapold = 0; 00425 int i; 00426 00427 for (i = 0; i < n; i++) { 00428 if (!cmap[i]) { 00429 cell[i] = 0; 00430 continue; 00431 } 00432 00433 if (cmap[i] == cmapold) { 00434 cell[i] = cell[i-1]; 00435 continue; 00436 } 00437 00438 cell[i] = data[cmap[i] - 1]; 00439 00440 cmapold = cmap[i]; 00441 } 00442 } 00443 00444 /*--------------------------------------------------------------------------*/ 00445 00446 static void gdal_values_double(int fd, const double *data, 00447 const COLUMN_MAPPING *cmap, int nbytes, 00448 DCELL *cell, int n) 00449 { 00450 COLUMN_MAPPING cmapold = 0; 00451 int i; 00452 00453 for (i = 0; i < n; i++) { 00454 if (!cmap[i]) { 00455 cell[i] = 0; 00456 continue; 00457 } 00458 00459 if (cmap[i] == cmapold) { 00460 cell[i] = cell[i-1]; 00461 continue; 00462 } 00463 00464 cell[i] = data[cmap[i] - 1]; 00465 00466 cmapold = cmap[i]; 00467 } 00468 } 00469 00470 /*--------------------------------------------------------------------------*/ 00471 00472 #endif 00473 00474 /*--------------------------------------------------------------------------*/ 00475 00476 /* transfer_to_cell_XY takes bytes from fcb->data, converts these bytes with 00477 the appropriate procedure (e.g. XDR or byte reordering) into type X 00478 values which are put into array G__.work_buf. 00479 finally the values in G__.work_buf are converted into 00480 type Y and put into 'cell'. 00481 if type X == type Y the intermediate step of storing the values in 00482 G__.work_buf might be ommited. check the appropriate function for XY to 00483 determine the procedure of conversion. 00484 */ 00485 00486 /*--------------------------------------------------------------------------*/ 00487 00488 static void transfer_to_cell_XX(int fd, void *cell) 00489 { 00490 static void (*cell_values_type[3]) () = { 00491 cell_values_int, cell_values_float, cell_values_double}; 00492 #ifdef HAVE_GDAL 00493 static void (*gdal_values_type[3]) () = { 00494 gdal_values_int, gdal_values_float, gdal_values_double}; 00495 #endif 00496 struct fileinfo *fcb = &G__.fileinfo[fd]; 00497 00498 #ifdef HAVE_GDAL 00499 if (fcb->gdal) 00500 (gdal_values_type[fcb->map_type]) (fd, fcb->data, fcb->col_map, 00501 fcb->cur_nbytes, cell, 00502 G__.window.cols); 00503 else 00504 #endif 00505 (cell_values_type[fcb->map_type]) (fd, fcb->data, fcb->col_map, 00506 fcb->cur_nbytes, cell, 00507 G__.window.cols); 00508 } 00509 00510 /*--------------------------------------------------------------------------*/ 00511 00512 static void transfer_to_cell_fi(int fd, void *cell) 00513 { 00514 struct fileinfo *fcb = &G__.fileinfo[fd]; 00515 int i; 00516 00517 transfer_to_cell_XX(fd, G__.work_buf); 00518 00519 for (i = 0; i < G__.window.cols; i++) 00520 ((CELL *) cell)[i] = (fcb->col_map[i] == 0) 00521 ? 0 00522 : G_quant_get_cell_value(&fcb->quant, 00523 ((FCELL *) G__.work_buf)[i]); 00524 } 00525 00526 static void transfer_to_cell_di(int fd, void *cell) 00527 { 00528 struct fileinfo *fcb = &G__.fileinfo[fd]; 00529 int i; 00530 00531 transfer_to_cell_XX(fd, G__.work_buf); 00532 00533 for (i = 0; i < G__.window.cols; i++) 00534 ((CELL *) cell)[i] = (fcb->col_map[i] == 0) 00535 ? 0 00536 : G_quant_get_cell_value(&fcb->quant, 00537 ((DCELL *) G__.work_buf)[i]); 00538 } 00539 00540 /*--------------------------------------------------------------------------*/ 00541 00542 static void transfer_to_cell_if(int fd, void *cell) 00543 { 00544 int i; 00545 00546 transfer_to_cell_XX(fd, G__.work_buf); 00547 00548 for (i = 0; i < G__.window.cols; i++) 00549 ((FCELL *) cell)[i] = ((CELL *) G__.work_buf)[i]; 00550 } 00551 00552 static void transfer_to_cell_df(int fd, void *cell) 00553 { 00554 int i; 00555 00556 transfer_to_cell_XX(fd, G__.work_buf); 00557 00558 for (i = 0; i < G__.window.cols; i++) 00559 ((FCELL *) cell)[i] = ((DCELL *) G__.work_buf)[i]; 00560 } 00561 00562 /*--------------------------------------------------------------------------*/ 00563 00564 static void transfer_to_cell_id(int fd, void *cell) 00565 { 00566 int i; 00567 00568 transfer_to_cell_XX(fd, G__.work_buf); 00569 00570 for (i = 0; i < G__.window.cols; i++) 00571 ((DCELL *) cell)[i] = ((CELL *) G__.work_buf)[i]; 00572 } 00573 00574 static void transfer_to_cell_fd(int fd, void *cell) 00575 { 00576 int i; 00577 00578 transfer_to_cell_XX(fd, G__.work_buf); 00579 00580 for (i = 0; i < G__.window.cols; i++) 00581 ((DCELL *) cell)[i] = ((FCELL *) G__.work_buf)[i]; 00582 } 00583 00584 /*--------------------------------------------------------------------------*/ 00585 00586 /*--------------------------------------------------------------------------*/ 00587 00588 /*--------------------------------------------------------------------------*/ 00589 /* 00590 * works for all map types and doesn't consider 00591 * null row corresponding to the requested row 00592 */ 00593 static int get_map_row_nomask(int fd, void *rast, int row, 00594 RASTER_MAP_TYPE data_type) 00595 { 00596 static void (*transfer_to_cell_FtypeOtype[3][3]) () = { { 00597 transfer_to_cell_XX, transfer_to_cell_if, transfer_to_cell_id}, { 00598 transfer_to_cell_fi, transfer_to_cell_XX, transfer_to_cell_fd}, { 00599 transfer_to_cell_di, transfer_to_cell_df, transfer_to_cell_XX}}; 00600 struct fileinfo *fcb = &G__.fileinfo[fd]; 00601 int r; 00602 int rowStatus; 00603 00604 rowStatus = compute_window_row(fd, row, &r); 00605 00606 if (rowStatus <= 0) { 00607 fcb->cur_row = -1; 00608 G_zero_raster_buf(rast, data_type); 00609 return rowStatus; 00610 } 00611 00612 /* read cell file row if not in memory */ 00613 if (r != fcb->cur_row) { 00614 fcb->cur_row = r; 00615 00616 if (read_data(fd, fcb->cur_row, fcb->data, &fcb->cur_nbytes) < 0) { 00617 G_zero_raster_buf(rast, data_type); 00618 00619 if (!fcb->io_error) { 00620 if (fcb->cellhd.compressed) 00621 G_warning(_("Error reading compressed map <%s@%s>, row %d"), 00622 fcb->name, fcb->mapset, r); 00623 else 00624 G_warning(_("Error reading map <%s@%s>, row %d"), 00625 fcb->name, fcb->mapset, r); 00626 00627 fcb->io_error = 1; 00628 } 00629 return -1; 00630 } 00631 } 00632 00633 (transfer_to_cell_FtypeOtype[fcb->map_type][data_type]) (fd, rast); 00634 00635 return 1; 00636 } 00637 00638 /*--------------------------------------------------------------------------*/ 00639 00640 static int get_map_row_no_reclass(int fd, void *rast, int row, 00641 RASTER_MAP_TYPE data_type, int null_is_zero, 00642 int with_mask) 00643 { 00644 int stat; 00645 00646 stat = get_map_row_nomask(fd, rast, row, data_type); 00647 if (stat < 0) 00648 return stat; 00649 00650 stat = embed_nulls(fd, rast, row, data_type, null_is_zero, with_mask); 00651 if (stat < 0) 00652 return stat; 00653 00654 return 1; 00655 } 00656 00657 /*--------------------------------------------------------------------------*/ 00658 00659 static int get_map_row(int fd, void *rast, int row, RASTER_MAP_TYPE data_type, 00660 int null_is_zero, int with_mask) 00661 { 00662 struct fileinfo *fcb = &G__.fileinfo[fd]; 00663 int size = G_raster_size(data_type); 00664 void *buf; 00665 int type; 00666 int stat; 00667 int i; 00668 00669 if (fcb->reclass_flag && data_type != CELL_TYPE) { 00670 buf = G__.temp_buf; 00671 type = CELL_TYPE; 00672 } 00673 else { 00674 buf = rast; 00675 type = data_type; 00676 } 00677 00678 stat = 00679 get_map_row_no_reclass(fd, buf, row, type, null_is_zero, with_mask); 00680 if (stat < 0) 00681 return stat; 00682 00683 if (!fcb->reclass_flag) 00684 return 1; 00685 00686 /* if the map is reclass table, get and 00687 reclass CELL row and copy results to needed type */ 00688 00689 do_reclass_int(fd, buf, null_is_zero); 00690 00691 if (data_type == CELL_TYPE) 00692 return 1; 00693 00694 for (i = 0; i < G__.window.cols; i++) { 00695 G_set_raster_value_c(rast, G__.temp_buf[i], data_type); 00696 rast = G_incr_void_ptr(rast, size); 00697 } 00698 00699 return 1; 00700 } 00701 00702 /*--------------------------------------------------------------------------*/ 00703 00704 /*--------------------------------------------------------------------------*/ 00705 00706 /*--------------------------------------------------------------------------*/ 00707 00740 int G_get_map_row_nomask(int fd, CELL * buf, int row) 00741 { 00742 return get_map_row(fd, buf, row, CELL_TYPE, 1, 0); 00743 } 00744 00760 int G_get_raster_row_nomask(int fd, void *buf, int row, 00761 RASTER_MAP_TYPE data_type) 00762 { 00763 return get_map_row(fd, buf, row, data_type, 0, 0); 00764 } 00765 00781 int G_get_c_raster_row_nomask(int fd, CELL * buf, int row) 00782 { 00783 return G_get_raster_row_nomask(fd, buf, row, CELL_TYPE); 00784 } 00785 00801 int G_get_f_raster_row_nomask(int fd, FCELL * buf, int row) 00802 { 00803 return G_get_raster_row_nomask(fd, buf, row, FCELL_TYPE); 00804 } 00805 00821 int G_get_d_raster_row_nomask(int fd, DCELL * buf, int row) 00822 { 00823 return G_get_raster_row_nomask(fd, buf, row, DCELL_TYPE); 00824 } 00825 00826 /*--------------------------------------------------------------------------*/ 00827 00848 int G_get_map_row(int fd, CELL * buf, int row) 00849 { 00850 return get_map_row(fd, buf, row, CELL_TYPE, 1, 1); 00851 } 00852 00896 int G_get_raster_row(int fd, void *buf, int row, RASTER_MAP_TYPE data_type) 00897 { 00898 return get_map_row(fd, buf, row, data_type, 0, 1); 00899 } 00900 00923 int G_get_c_raster_row(int fd, CELL * buf, int row) 00924 { 00925 return G_get_raster_row(fd, buf, row, CELL_TYPE); 00926 } 00927 00946 int G_get_f_raster_row(int fd, FCELL * buf, int row) 00947 { 00948 return G_get_raster_row(fd, buf, row, FCELL_TYPE); 00949 } 00950 00966 int G_get_d_raster_row(int fd, DCELL * buf, int row) 00967 { 00968 return G_get_raster_row(fd, buf, row, DCELL_TYPE); 00969 } 00970 00971 /*--------------------------------------------------------------------------*/ 00972 00973 /*--------------------------------------------------------------------------*/ 00974 00975 /*--------------------------------------------------------------------------*/ 00976 00977 static int open_null_read(int fd) 00978 { 00979 struct fileinfo *fcb = &G__.fileinfo[fd]; 00980 char *name, *mapset, *dummy; 00981 int null_fd; 00982 00983 if (fcb->null_file_exists == 0) 00984 return -1; 00985 00986 if (fcb->reclass_flag) { 00987 name = fcb->reclass.name; 00988 mapset = fcb->reclass.mapset; 00989 } 00990 else { 00991 name = fcb->name; 00992 mapset = fcb->mapset; 00993 } 00994 00995 dummy = G_find_file_misc("cell_misc", NULL_FILE, name, mapset); 00996 00997 if (!dummy) { 00998 /* G_warning("unable to find [%s]",path); */ 00999 fcb->null_file_exists = 0; 01000 return -1; 01001 } 01002 01003 null_fd = G_open_old_misc("cell_misc", NULL_FILE, name, mapset); 01004 if (null_fd < 0) 01005 return -1; 01006 01007 fcb->null_file_exists = 1; 01008 01009 return null_fd; 01010 } 01011 01012 static int read_null_bits(int null_fd, unsigned char *flags, int row, 01013 int cols, int fd) 01014 { 01015 off_t offset; 01016 ssize_t size; 01017 int R; 01018 01019 if (compute_window_row(fd, row, &R) <= 0) { 01020 G__init_null_bits(flags, cols); 01021 return 1; 01022 } 01023 01024 if (null_fd < 0) 01025 return -1; 01026 01027 size = G__null_bitstream_size(cols); 01028 offset = (off_t) size *R; 01029 01030 if (lseek(null_fd, offset, SEEK_SET) < 0) { 01031 G_warning(_("Error reading null row %d"), R); 01032 return -1; 01033 } 01034 01035 if (read(null_fd, flags, size) != size) { 01036 G_warning(_("Error reading null row %d"), R); 01037 return -1; 01038 } 01039 01040 return 1; 01041 } 01042 01043 static void get_null_value_row_nomask(int fd, char *flags, int row) 01044 { 01045 struct fileinfo *fcb = &G__.fileinfo[fd]; 01046 int i, j, null_fd; 01047 01048 if (row > G__.window.rows || row < 0) { 01049 G_warning(_("Reading raster map <%s@%s> request for row %d is outside region"), 01050 fcb->name, fcb->mapset, row); 01051 } 01052 01053 if ((fcb->min_null_row > row) || 01054 (fcb->min_null_row + NULL_ROWS_INMEM - 1 < row)) 01055 /* the null row row is not in memory */ 01056 { 01057 /* read in NULL_ROWS_INMEM rows from null file 01058 so that the requested row is between fcb->min_null_row 01059 and fcb->min_null_row + NULL_ROWS_INMEM */ 01060 01061 fcb->min_null_row = (row / NULL_ROWS_INMEM) * NULL_ROWS_INMEM; 01062 01063 null_fd = open_null_read(fd); 01064 01065 for (i = 0; i < NULL_ROWS_INMEM; i++) { 01066 /* G__.window.rows doesn't have to be a multiple of NULL_ROWS_INMEM */ 01067 if (i + fcb->min_null_row >= G__.window.rows) 01068 break; 01069 01070 if (read_null_bits(null_fd, fcb->null_work_buf, 01071 i + fcb->min_null_row, fcb->cellhd.cols, 01072 fd) < 0) { 01073 if (fcb->map_type == CELL_TYPE) { 01074 /* 01075 If can't read null row, assume that all map 0's are nulls 01076 use allocated G__.mask_buf to read map row */ 01077 get_map_row_nomask(fd, (void *)G__.mask_buf, 01078 i + fcb->min_null_row, CELL_TYPE); 01079 for (j = 0; j < G__.window.cols; j++) { 01080 if (G__.mask_buf[j] == 0) 01081 flags[j] = 1; 01082 else 01083 flags[j] = 0; 01084 } 01085 } 01086 else { /* fp map */ 01087 01088 /* if can't read null row, assume that all data is valid */ 01089 G_zero(flags, sizeof(char) * G__.window.cols); 01090 /* the flags row is ready now */ 01091 } 01092 } /*if no null file */ 01093 else { 01094 /* copy null row to flags row translated by window column mapping */ 01095 /* the fcb->NULL_ROWS[row-fcb->min_null_row] has G__.window.cols bits, */ 01096 /* the fcb->null_work_buf has size fcb->cellhd.cols */ 01097 for (j = 0; j < G__.window.cols; j++) { 01098 if (!fcb->col_map[j]) 01099 flags[j] = 1; 01100 else 01101 flags[j] = G__check_null_bit(fcb->null_work_buf, 01102 fcb->col_map[j] - 1, 01103 fcb->cellhd.cols); 01104 } 01105 } 01106 /* remember the null row for i for the future reference */ 01107 01108 /*bf-We should take of the size - or we get 01109 zeros running on their own after flags convertions -A.Sh. */ 01110 fcb->NULL_ROWS[i] = G_realloc(fcb->NULL_ROWS[i], 01111 G__null_bitstream_size(G__.window. 01112 cols) + 1); 01113 if (fcb->NULL_ROWS[i] == NULL) 01114 G_fatal_error("get_null_value_row_nomask: %s", 01115 _("Unable to realloc buffer")); 01116 01117 G__convert_01_flags(flags, fcb->NULL_ROWS[i], G__.window.cols); 01118 01119 } /* for loop */ 01120 01121 if (null_fd > 0) 01122 close(null_fd); 01123 } /* row is not in memory */ 01124 01125 /* copy null file data translated by column mapping to user null row */ 01126 /* the user requested flags row is of size G__.window.cols */ 01127 G__convert_flags_01(flags, fcb->NULL_ROWS[row - fcb->min_null_row], 01128 G__.window.cols); 01129 } 01130 01131 /*--------------------------------------------------------------------------*/ 01132 01133 #ifdef HAVE_GDAL 01134 01135 static void get_null_value_row_gdal(int fd, char *flags, int row) 01136 { 01137 struct fileinfo *fcb = &G__.fileinfo[fd]; 01138 DCELL *tmp_buf = G_allocate_d_raster_buf(); 01139 int i; 01140 01141 if (get_map_row_nomask(fd, tmp_buf, row, DCELL_TYPE) <= 0) { 01142 memset(flags, 1, G__.window.cols); 01143 G_free(tmp_buf); 01144 return; 01145 } 01146 01147 for (i = 0; i < G__.window.cols; i++) 01148 /* note: using == won't work if the null value is NaN */ 01149 flags[i] = !fcb->col_map[i] || 01150 memcmp(&tmp_buf[i], &fcb->gdal->null_val, sizeof(DCELL)) == 0; 01151 01152 G_free(tmp_buf); 01153 } 01154 01155 #endif 01156 01157 /*--------------------------------------------------------------------------*/ 01158 01159 /*--------------------------------------------------------------------------*/ 01160 01161 static void embed_mask(char *flags, int row) 01162 { 01163 int i; 01164 01165 if (G__.auto_mask <= 0) 01166 return; 01167 01168 if (get_map_row_nomask(G__.mask_fd, G__.mask_buf, row, CELL_TYPE) < 0) 01169 return; 01170 01171 if (G__.fileinfo[G__.mask_fd].reclass_flag) 01172 do_reclass_int(G__.mask_fd, G__.mask_buf, 1); 01173 01174 for (i = 0; i < G__.window.cols; i++) 01175 if (G__.mask_buf[i] == 0) 01176 flags[i] = 1; 01177 } 01178 01179 static void get_null_value_row(int fd, char *flags, int row, int with_mask) 01180 { 01181 #ifdef HAVE_GDAL 01182 struct fileinfo *fcb = &G__.fileinfo[fd]; 01183 if (fcb->gdal) 01184 get_null_value_row_gdal(fd, flags, row); 01185 else 01186 #endif 01187 get_null_value_row_nomask(fd, flags, row); 01188 01189 if (with_mask) 01190 embed_mask(flags, row); 01191 } 01192 01193 static int embed_nulls(int fd, void *buf, int row, RASTER_MAP_TYPE map_type, 01194 int null_is_zero, int with_mask) 01195 { 01196 struct fileinfo *fcb = &G__.fileinfo[fd]; 01197 int i; 01198 01199 /* this is because without null file the nulls can be only due to 0's 01200 in data row or mask */ 01201 if (null_is_zero && !fcb->null_file_exists 01202 && (G__.auto_mask <= 0 || !with_mask)) 01203 return 1; 01204 01205 get_null_value_row(fd, G__.null_buf, row, with_mask); 01206 01207 for (i = 0; i < G__.window.cols; i++) { 01208 /* also check for nulls which might be already embedded by quant 01209 rules in case of fp map. */ 01210 if (G__.null_buf[i] || G_is_null_value(buf, map_type)) { 01211 /* G__set_[f/d]_null_value() sets it to 0 is the embedded mode 01212 is not set and calls G_set_[f/d]_null_value() otherwise */ 01213 G__set_null_value(buf, 1, null_is_zero, map_type); 01214 } 01215 buf = G_incr_void_ptr(buf, G_raster_size(map_type)); 01216 } 01217 01218 return 1; 01219 } 01220 01221 /*--------------------------------------------------------------------------*/ 01222 01223 /*--------------------------------------------------------------------------*/ 01224 01225 /*--------------------------------------------------------------------------*/ 01226 01243 int G_get_null_value_row(int fd, char *flags, int row) 01244 { 01245 get_null_value_row(fd, flags, row, 1); 01246 01247 return 1; 01248 }