~ubuntu-branches/ubuntu/wily/grass/wily

« back to all changes in this revision

Viewing changes to lib/raster/range.c

Tags: 7.0.0~rc1+ds1-1~exp1
* New upstream release candidate.
* Repack upstream tarball, remove precompiled Python objects.
* Add upstream metadata.
* Update gbp.conf and Vcs-Git URL to use the experimental branch.
* Update watch file for GRASS 7.0.
* Drop build dependencies for Tcl/Tk, add build dependencies:
  python-numpy, libnetcdf-dev, netcdf-bin, libblas-dev, liblapack-dev
* Update Vcs-Browser URL to use cgit instead of gitweb.
* Update paths to use grass70.
* Add configure options: --with-netcdf, --with-blas, --with-lapack,
  remove --with-tcltk-includes.
* Update patches for GRASS 7.
* Update copyright file, changes:
  - Update copyright years
  - Group files by license
  - Remove unused license sections
* Add patches for various typos.
* Fix desktop file with patch instead of d/rules.
* Use minimal dh rules.
* Bump Standards-Version to 3.9.6, no changes.
* Use dpkg-maintscript-helper to replace directories with symlinks.
  (closes: #776349)
* Update my email to use @debian.org address.

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/*!
 
2
 * \file lib/raster/range.c
 
3
 *
 
4
 * \brief Raster Library - Raster range file management
 
5
 *
 
6
 * (C) 2001-2009 GRASS Development Team
 
7
 *
 
8
 * This program is free software under the GNU General Public License 
 
9
 * (>=v2). Read the file COPYING that comes with GRASS for details.
 
10
 *
 
11
 * \author Original author CERL
 
12
 */
 
13
 
 
14
#include <unistd.h>
 
15
#include <rpc/types.h>          /* need this for sgi */
 
16
 
 
17
#include <grass/raster.h>
 
18
#include <grass/glocale.h>
 
19
 
 
20
#include "R.h"
 
21
 
 
22
#define DEFAULT_CELL_MIN 1
 
23
#define DEFAULT_CELL_MAX 255
 
24
 
 
25
/*!
 
26
   \brief Remove floating-point range
 
27
 
 
28
   Note: For internal use only.
 
29
 
 
30
   \param name map name
 
31
 */
 
32
void Rast__remove_fp_range(const char *name)
 
33
{
 
34
    G_remove_misc("cell_misc", "f_range", name);
 
35
}
 
36
 
 
37
/*!
 
38
 * \brief Construct default range
 
39
 *
 
40
 * Sets the integer range to [1,255]
 
41
 *
 
42
 * \param[out] r pointer to Range structure which holds range info
 
43
 */
 
44
void Rast_construct_default_range(struct Range *range)
 
45
{
 
46
    Rast_update_range(DEFAULT_CELL_MIN, range);
 
47
    Rast_update_range(DEFAULT_CELL_MAX, range);
 
48
}
 
49
 
 
50
/*!
 
51
 * \brief Read floating-point range
 
52
 *
 
53
 * Read the floating point range file <i>drange</i>. This file is
 
54
 * written in binary using XDR format.
 
55
 *
 
56
 * An empty range file indicates that the min, max are undefined. This
 
57
 * is a valid case, and the result should be an initialized range
 
58
 * struct with no defined min/max.  If the range file is missing and
 
59
 * the map is a floating-point map, this function will create a
 
60
 * default range by calling G_construct_default_range().
 
61
 *
 
62
 * \param name map name
 
63
 * \param mapset mapset name
 
64
 * \param drange pointer to FPRange structure which holds fp range
 
65
 *
 
66
 * \return 1 on success
 
67
 * \return 2 range is empty
 
68
 * \return -1 on error
 
69
 */
 
70
int Rast_read_fp_range(const char *name, const char *mapset,
 
71
                       struct FPRange *drange)
 
72
{
 
73
    struct Range range;
 
74
    int fd;
 
75
    char xdr_buf[2][XDR_DOUBLE_NBYTES];
 
76
    DCELL dcell1, dcell2;
 
77
 
 
78
    Rast_init();
 
79
    Rast_init_fp_range(drange);
 
80
 
 
81
    if (Rast_map_type(name, mapset) == CELL_TYPE) {
 
82
        /* if map is integer
 
83
           read integer range and convert it to double */
 
84
 
 
85
        if (Rast_read_range(name, mapset, &range) >= 0) {
 
86
            /* if the integer range is empty */
 
87
            if (range.first_time)
 
88
                return 2;
 
89
 
 
90
            Rast_update_fp_range((DCELL) range.min, drange);
 
91
            Rast_update_fp_range((DCELL) range.max, drange);
 
92
            return 1;
 
93
        }
 
94
        return -1;
 
95
    }
 
96
 
 
97
    fd = -1;
 
98
 
 
99
    if (G_find_file2_misc("cell_misc", "f_range", name, mapset)) {
 
100
        fd = G_open_old_misc("cell_misc", "f_range", name, mapset);
 
101
        if (fd < 0) {
 
102
            G_warning(_("Unable to read fp range file for <%s>"),
 
103
                      G_fully_qualified_name(name, mapset));
 
104
            return -1;
 
105
        }
 
106
 
 
107
        if (read(fd, xdr_buf, sizeof(xdr_buf)) != sizeof(xdr_buf)) {
 
108
            /* if the f_range file exists, but empty file, meaning Nulls */
 
109
            close(fd);
 
110
            G_debug(1, "Empty fp range file meaning Nulls for <%s>",
 
111
                      G_fully_qualified_name(name, mapset));
 
112
            return 2;
 
113
        }
 
114
 
 
115
        G_xdr_get_double(&dcell1, xdr_buf[0]);
 
116
        G_xdr_get_double(&dcell2, xdr_buf[1]);
 
117
 
 
118
        Rast_update_fp_range(dcell1, drange);
 
119
        Rast_update_fp_range(dcell2, drange);
 
120
        close(fd);
 
121
    }
 
122
 
 
123
    return 1;
 
124
}
 
125
 
 
126
/*!
 
127
 * \brief Read raster range (CELL)
 
128
 *
 
129
 * This routine reads the range information for the raster map
 
130
 * <i>name</i> in <i>mapset</i> into the <i>range</i> structure.
 
131
 *
 
132
 * A diagnostic message is printed and -1 is returned if there is an error
 
133
 * reading the range file. Otherwise, 0 is returned.
 
134
 *
 
135
 * Old range file (those with 4 numbers) should treat zeros in this
 
136
 * file as NULL-values. New range files (those with just 2 numbers)
 
137
 * should treat these numbers as real data (zeros are real data in
 
138
 * this case).  An empty range file indicates that the min, max are
 
139
 * undefined. This is a valid case, and the result should be an
 
140
 * initialized range struct with no defined min/max. If the range file
 
141
 * is missing and the map is a floating-point map, this function will
 
142
 * create a default range by calling G_construct_default_range().
 
143
 *
 
144
 * \param name map name
 
145
 * \param mapset mapset name
 
146
 * \param[out] range pointer to Range structure which holds range info
 
147
 *
 
148
 * \return -1 on error
 
149
 * \return 1 on success
 
150
 * \return 2 if range is empty
 
151
 * \return 3 if raster map is floating-point, get range from quant rules
 
152
 */
 
153
int Rast_read_range(const char *name, const char *mapset, struct Range *range)
 
154
{
 
155
    FILE *fd;
 
156
    CELL x[4];
 
157
    char buf[200];
 
158
    int n, count;
 
159
    struct Quant quant;
 
160
    struct FPRange drange;
 
161
    Rast_init_range(range);
 
162
    fd = NULL;
 
163
 
 
164
    /* if map is not integer, read quant rules, and get limits */
 
165
    if (Rast_map_type(name, mapset) != CELL_TYPE) {
 
166
        DCELL dmin, dmax;
 
167
 
 
168
        if (Rast_read_quant(name, mapset, &quant) < 0) {
 
169
            G_warning(_("Unable to read quant rules for raster map <%s>"),
 
170
                      G_fully_qualified_name(name, mapset));
 
171
            return -1;
 
172
        }
 
173
        if (Rast_quant_is_truncate(&quant) || Rast_quant_is_round(&quant)) {
 
174
            if (Rast_read_fp_range(name, mapset, &drange) >= 0) {
 
175
                Rast_get_fp_range_min_max(&drange, &dmin, &dmax);
 
176
                if (Rast_quant_is_truncate(&quant)) {
 
177
                    x[0] = (CELL) dmin;
 
178
                    x[1] = (CELL) dmax;
 
179
                }
 
180
                else {          /* round */
 
181
 
 
182
                    if (dmin > 0)
 
183
                        x[0] = (CELL) (dmin + .5);
 
184
                    else
 
185
                        x[0] = (CELL) (dmin - .5);
 
186
                    if (dmax > 0)
 
187
                        x[1] = (CELL) (dmax + .5);
 
188
                    else
 
189
                        x[1] = (CELL) (dmax - .5);
 
190
                }
 
191
            }
 
192
            else
 
193
                return -1;
 
194
        }
 
195
        else
 
196
            Rast_quant_get_limits(&quant, &dmin, &dmax, &x[0], &x[1]);
 
197
 
 
198
        Rast_update_range(x[0], range);
 
199
        Rast_update_range(x[1], range);
 
200
        return 3;
 
201
    }
 
202
 
 
203
    if (G_find_file2_misc("cell_misc", "range", name, mapset)) {
 
204
        fd = G_fopen_old_misc("cell_misc", "range", name, mapset);
 
205
        if (!fd) {
 
206
            G_warning(_("Unable to read range file for <%s>"),
 
207
                      G_fully_qualified_name(name, mapset));
 
208
            return -1;
 
209
        }
 
210
 
 
211
        /* if range file exists but empty */
 
212
        if (!fgets(buf, sizeof buf, fd)) {
 
213
            if (fd)
 
214
                fclose(fd);
 
215
            return 2;
 
216
        }
 
217
 
 
218
        x[0] = x[1] = x[2] = x[3] = 0;
 
219
        count = sscanf(buf, "%d%d%d%d", &x[0], &x[1], &x[2], &x[3]);
 
220
 
 
221
        /* if wrong format */
 
222
        if (count <= 0) {
 
223
            if (fd)
 
224
                fclose(fd);
 
225
 
 
226
            G_warning(_("Unable to read range file for <%s>"),
 
227
                      G_fully_qualified_name(name, mapset));
 
228
            return -1;
 
229
        }
 
230
 
 
231
        for (n = 0; n < count; n++) {
 
232
            /* if count==4, the range file is old (4.1) and 0's in it
 
233
               have to be ignored */
 
234
            if (count < 4 || x[n])
 
235
                Rast_update_range((CELL) x[n], range);
 
236
        }
 
237
        fclose(fd);
 
238
    }
 
239
 
 
240
    return 1;
 
241
}
 
242
 
 
243
/*!
 
244
 * \brief Write raster range file
 
245
 *
 
246
 * This routine writes the range information for the raster map
 
247
 * <i>name</i> in the current mapset from the <i>range</i> structure.
 
248
 * A diagnostic message is printed and -1 is returned if there is an
 
249
 * error writing the range file. Otherwise, 0 is returned.
 
250
 *
 
251
 * This routine only writes 2 numbers (min,max) to the range
 
252
 * file, instead of the 4 (pmin,pmax,nmin,nmax) previously written.
 
253
 * If there is no defined min,max, an empty file is written.
 
254
 *
 
255
 * \param name map name
 
256
 * \param range pointer to Range structure which holds range info
 
257
 */
 
258
void Rast_write_range(const char *name, const struct Range *range)
 
259
{
 
260
    FILE *fp;
 
261
 
 
262
    if (Rast_map_type(name, G_mapset()) != CELL_TYPE) {
 
263
        G_remove_misc("cell_misc", "range", name);      /* remove the old file with this name */
 
264
        G_fatal_error(_("Unable to write range file for <%s>"), name);
 
265
    }
 
266
 
 
267
    fp = G_fopen_new_misc("cell_misc", "range", name);
 
268
    if (!fp) {
 
269
        G_remove_misc("cell_misc", "range", name);      /* remove the old file with this name */
 
270
        G_fatal_error(_("Unable to write range file for <%s>"), name);
 
271
    }
 
272
 
 
273
    /* if range has been updated */
 
274
    if (!range->first_time)
 
275
        fprintf(fp, "%ld %ld\n", (long)range->min, (long)range->max);
 
276
 
 
277
    fclose(fp);
 
278
}
 
279
 
 
280
/*!
 
281
 * \brief Write raster range file (floating-point)
 
282
 *
 
283
 * Write the floating point range file <tt>f_range</tt>. This file is
 
284
 * written in binary using XDR format. If there is no defined min/max
 
285
 * in <em>range</em>, an empty <tt>f_range</tt> file is created.
 
286
 *
 
287
 * \param name map name
 
288
 * \param range pointer to FPRange which holds fp range info
 
289
 */
 
290
void Rast_write_fp_range(const char *name, const struct FPRange *range)
 
291
{
 
292
    int fd;
 
293
    char xdr_buf[2][XDR_DOUBLE_NBYTES];
 
294
 
 
295
    Rast_init();
 
296
 
 
297
    fd = G_open_new_misc("cell_misc", "f_range", name);
 
298
    if (fd < 0) {
 
299
        G_remove_misc("cell_misc", "f_range", name);
 
300
        G_fatal_error(_("Unable to write range file for <%s>"), name);
 
301
    }
 
302
 
 
303
    /* if range hasn't been updated, write empty file meaning Nulls */
 
304
    if (range->first_time) {
 
305
        close(fd);
 
306
        return;
 
307
    }
 
308
 
 
309
    G_xdr_put_double(xdr_buf[0], &range->min);
 
310
    G_xdr_put_double(xdr_buf[1], &range->max);
 
311
 
 
312
    if (write(fd, xdr_buf, sizeof(xdr_buf)) != sizeof(xdr_buf)) {
 
313
        G_remove_misc("cell_misc", "f_range", name);
 
314
        G_fatal_error(_("Unable to write range file for <%s>"), name);
 
315
    }
 
316
 
 
317
    close(fd);
 
318
}
 
319
 
 
320
/*!
 
321
 * \brief Update range structure (CELL)
 
322
 *
 
323
 * Compares the <i>cat</i> value with the minimum and maximum values
 
324
 * in the <i>range</i> structure, modifying the range if <i>cat</i>
 
325
 * extends the range.
 
326
 *
 
327
 * NULL-values must be detected and ignored.
 
328
 *
 
329
 * \param cat raster value
 
330
 * \param range pointer to Range structure which holds range info
 
331
 */
 
332
void Rast_update_range(CELL cat, struct Range *range)
 
333
{
 
334
    if (!Rast_is_c_null_value(&cat)) {
 
335
        if (range->first_time) {
 
336
            range->first_time = 0;
 
337
            range->min = cat;
 
338
            range->max = cat;
 
339
            return;
 
340
        }
 
341
        if (cat < range->min)
 
342
            range->min = cat;
 
343
        if (cat > range->max)
 
344
            range->max = cat;
 
345
    }
 
346
}
 
347
 
 
348
/*!
 
349
 * \brief Update range structure (floating-point)
 
350
 *
 
351
 * Compares the <i>cat</i> value with the minimum and maximum values
 
352
 * in the <i>range</i> structure, modifying the range if <i>cat</i>
 
353
 * extends the range.
 
354
 *
 
355
 * NULL-values must be detected and ignored.
 
356
 *
 
357
 * \param val raster value
 
358
 * \param range pointer to Range structure which holds range info
 
359
 */
 
360
void Rast_update_fp_range(DCELL val, struct FPRange *range)
 
361
{
 
362
    if (!Rast_is_d_null_value(&val)) {
 
363
        if (range->first_time) {
 
364
            range->first_time = 0;
 
365
            range->min = val;
 
366
            range->max = val;
 
367
            return;
 
368
        }
 
369
        if (val < range->min)
 
370
            range->min = val;
 
371
        if (val > range->max)
 
372
            range->max = val;
 
373
    }
 
374
}
 
375
 
 
376
/*!
 
377
 * \brief Update range structure based on raster row (CELL)
 
378
 *
 
379
 * This routine updates the <i>range</i> data just like
 
380
 * Rast_update_range(), but for <i>n</i> values from the <i>cell</i>
 
381
 * array.
 
382
 *
 
383
 * \param cell raster values
 
384
 * \param n number of values
 
385
 * \param range pointer to Range structure which holds range info
 
386
 */
 
387
void Rast_row_update_range(const CELL * cell, int n, struct Range *range)
 
388
{
 
389
    Rast__row_update_range(cell, n, range, 0);
 
390
}
 
391
 
 
392
/*!
 
393
 * \brief Update range structure based on raster row
 
394
 *
 
395
 * Note: for internal use only.
 
396
 *
 
397
 * \param cell raster values
 
398
 * \param n number of values
 
399
 * \param range pointer to Range structure which holds range info
 
400
 * \param ignore_zeros ignore zeros
 
401
 */
 
402
void Rast__row_update_range(const CELL * cell, int n,
 
403
                            struct Range *range, int ignore_zeros)
 
404
{
 
405
    CELL cat;
 
406
 
 
407
    while (n-- > 0) {
 
408
        cat = *cell++;
 
409
        if (Rast_is_c_null_value(&cat) || (ignore_zeros && !cat))
 
410
            continue;
 
411
        if (range->first_time) {
 
412
            range->first_time = 0;
 
413
            range->min = cat;
 
414
            range->max = cat;
 
415
            continue;
 
416
        }
 
417
        if (cat < range->min)
 
418
            range->min = cat;
 
419
        if (cat > range->max)
 
420
            range->max = cat;
 
421
    }
 
422
}
 
423
 
 
424
/*!
 
425
 * \brief Update range structure based on raster row (floating-point)
 
426
 *
 
427
 * This routine updates the <i>range</i> data just like
 
428
 * Rast_update_range(), but for <i>n</i> values from the <i>cell</i>
 
429
 * array.
 
430
 *
 
431
 * \param cell raster values
 
432
 * \param n number of values
 
433
 * \param range pointer to Range structure which holds range info
 
434
 * \param data_type raster type (CELL, FCELL, DCELL)
 
435
 */
 
436
void Rast_row_update_fp_range(const void *rast, int n,
 
437
                              struct FPRange *range,
 
438
                              RASTER_MAP_TYPE data_type)
 
439
{
 
440
    size_t size = Rast_cell_size(data_type);
 
441
    DCELL val = 0.0;
 
442
 
 
443
    while (n-- > 0) {
 
444
        switch (data_type) {
 
445
        case CELL_TYPE:
 
446
            val = (DCELL) * ((CELL *) rast);
 
447
            break;
 
448
        case FCELL_TYPE:
 
449
            val = (DCELL) * ((FCELL *) rast);
 
450
            break;
 
451
        case DCELL_TYPE:
 
452
            val = *((DCELL *) rast);
 
453
            break;
 
454
        }
 
455
 
 
456
        if (Rast_is_null_value(rast, data_type)) {
 
457
            rast = G_incr_void_ptr(rast, size);
 
458
            continue;
 
459
        }
 
460
        if (range->first_time) {
 
461
            range->first_time = 0;
 
462
            range->min = val;
 
463
            range->max = val;
 
464
        }
 
465
        else {
 
466
            if (val < range->min)
 
467
                range->min = val;
 
468
            if (val > range->max)
 
469
                range->max = val;
 
470
        }
 
471
 
 
472
        rast = G_incr_void_ptr(rast, size);
 
473
    }
 
474
}
 
475
 
 
476
/*!
 
477
 * \brief Initialize range structure
 
478
 *
 
479
 * Initializes the <i>range</i> structure for updates by
 
480
 * Rast_update_range() and Rast_row_update_range().
 
481
 *
 
482
 * Must set a flag in the range structure that indicates that no
 
483
 * min/max have been defined - probably a <tt>"first"</tt> boolean
 
484
 * flag.
 
485
 *
 
486
 * \param range pointer to Range structure which holds range info
 
487
 */
 
488
void Rast_init_range(struct Range *range)
 
489
{
 
490
    Rast_set_c_null_value(&(range->min), 1);
 
491
    Rast_set_c_null_value(&(range->max), 1);
 
492
    range->first_time = 1;
 
493
}
 
494
 
 
495
/*!
 
496
 * \brief Get range min and max
 
497
 *
 
498
 * The mininum and maximum CELL values are extracted from the
 
499
 * <i>range</i> structure.
 
500
 *
 
501
 * If the range structure has no defined min/max (first!=0) there will
 
502
 * not be a valid range. In this case the min and max returned must be
 
503
 * the NULL-value.
 
504
 *
 
505
 * \param range pointer to Range structure which holds range info
 
506
 * \param[out] min minimum value
 
507
 * \param[out] max maximum value
 
508
 */
 
509
void Rast_get_range_min_max(const struct Range *range, CELL * min, CELL * max)
 
510
{
 
511
    if (range->first_time) {
 
512
        Rast_set_c_null_value(min, 1);
 
513
        Rast_set_c_null_value(max, 1);
 
514
    }
 
515
    else {
 
516
        if (Rast_is_c_null_value(&(range->min)))
 
517
            Rast_set_c_null_value(min, 1);
 
518
        else
 
519
            *min = range->min;
 
520
 
 
521
        if (Rast_is_c_null_value(&(range->max)))
 
522
            Rast_set_c_null_value(max, 1);
 
523
        else
 
524
            *max = range->max;
 
525
    }
 
526
}
 
527
 
 
528
/*!
 
529
 * \brief Initialize fp range
 
530
 *
 
531
 * Must set a flag in the range structure that indicates that no
 
532
 * min/max have been defined - probably a <tt>"first"</tt> boolean
 
533
 * flag.
 
534
 *
 
535
 * \param range pointer to FPRange which holds fp range info
 
536
 */
 
537
void Rast_init_fp_range(struct FPRange *range)
 
538
{
 
539
    Rast_set_d_null_value(&(range->min), 1);
 
540
    Rast_set_d_null_value(&(range->max), 1);
 
541
    range->first_time = 1;
 
542
}
 
543
 
 
544
/*!
 
545
 * \brief Get minumum and maximum value from fp range
 
546
 *
 
547
 * Extract the min/max from the range structure <i>range</i>.  If the
 
548
 * range structure has no defined min/max (first!=0) there will not be
 
549
 * a valid range. In this case the min and max returned must be the
 
550
 * NULL-value.
 
551
 *
 
552
 * \param range pointer to FPRange which holds fp range info
 
553
 * \param[out] min minimum value
 
554
 * \param[out] max maximum value
 
555
 */
 
556
void Rast_get_fp_range_min_max(const struct FPRange *range,
 
557
                               DCELL * min, DCELL * max)
 
558
{
 
559
    if (range->first_time) {
 
560
        Rast_set_d_null_value(min, 1);
 
561
        Rast_set_d_null_value(max, 1);
 
562
    }
 
563
    else {
 
564
        if (Rast_is_d_null_value(&(range->min)))
 
565
            Rast_set_d_null_value(min, 1);
 
566
        else
 
567
            *min = range->min;
 
568
 
 
569
        if (Rast_is_d_null_value(&(range->max)))
 
570
            Rast_set_d_null_value(max, 1);
 
571
        else
 
572
            *max = range->max;
 
573
    }
 
574
}