[GRASS-SVN] r43433 - grass-addons/imagery/i.landsat.acca
svn_grass at osgeo.org
svn_grass at osgeo.org
Thu Sep 9 03:51:44 EDT 2010
Author: hamish
Date: 2010-09-09 07:51:43 +0000 (Thu, 09 Sep 2010)
New Revision: 43433
Modified:
grass-addons/imagery/i.landsat.acca/algorithm.c
grass-addons/imagery/i.landsat.acca/local_proto.h
grass-addons/imagery/i.landsat.acca/main.c
grass-addons/imagery/i.landsat.acca/tools.c
Log:
run tools/grass_indent.sh
Modified: grass-addons/imagery/i.landsat.acca/algorithm.c
===================================================================
--- grass-addons/imagery/i.landsat.acca/algorithm.c 2010-09-09 07:50:19 UTC (rev 43432)
+++ grass-addons/imagery/i.landsat.acca/algorithm.c 2010-09-09 07:51:43 UTC (rev 43433)
@@ -33,13 +33,13 @@
#define VARI 4
/*
- Con Landsat-7 funciona bien pero con Landsat-5 falla
- la primera pasada metiendo terreno desertico como nube fría
- y luego se equivoca en resolver los ambiguos en áreas
- con mucha humedad ambiental.
+ Con Landsat-7 funciona bien pero con Landsat-5 falla
+ la primera pasada metiendo terreno desertico como nube fría
+ y luego se equivoca en resolver los ambiguos en áreas
+ con mucha humedad ambiental.
- Habría que reajustar las constantes para Landsat-5
-*/
+ Habría que reajustar las constantes para Landsat-5
+ */
@@ -56,42 +56,40 @@
como opciones desde el programa main.
---------------------------------------------------------*/
-double th_1 = 0.08; /* Band 3 Brightness Threshold */
+double th_1 = 0.08; /* Band 3 Brightness Threshold */
double th_1_b = 0.07;
-double th_2[] = { -0.25, 0.70 }; /* Normalized Snow Difference Index */
+double th_2[] = { -0.25, 0.70 }; /* Normalized Snow Difference Index */
double th_2_b = 0.8;
-double th_3 = 300.; /* Band 6 Temperature Threshold */
-double th_4 = 225.; /* Band 5/6 Composite */
+double th_3 = 300.; /* Band 6 Temperature Threshold */
+double th_4 = 225.; /* Band 5/6 Composite */
double th_4_b = 0.08;
-double th_5 = 2.35; /* Band 4/3 Ratio */
-double th_6 = 2.16248; /* Band 4/2 Ratio */
-double th_7 = 1.0; /* Band 4/5 Ratio */;
-double th_8 = 210.; /* Band 5/6 Composite */
+double th_5 = 2.35; /* Band 4/3 Ratio */
+double th_6 = 2.16248; /* Band 4/2 Ratio */
+double th_7 = 1.0; /* Band 4/5 Ratio */ ;
+double th_8 = 210.; /* Band 5/6 Composite */
extern int hist_n;
#define K_BASE 240.
void acca_algorithm(int verbose, Gfile * out, Gfile band[],
- int two_pass, int with_shadow)
+ int two_pass, int with_shadow)
{
int i, count[5], hist_cold[hist_n], hist_warm[hist_n];
double x, value[5], signa[5], idesert, warm_ambiguous, shift;
/* Initialized varibles ... */
- for (i = 0; i < 5; i++)
- {
- count[i] = 0;
- value[i] = 0.;
+ for (i = 0; i < 5; i++) {
+ count[i] = 0;
+ value[i] = 0.;
}
- for (i = 0; i < hist_n; i++)
- {
- hist_cold[i] = hist_warm[i] = 0;
+ for (i = 0; i < hist_n; i++) {
+ hist_cold[i] = hist_warm[i] = 0;
}
/* FIRST FILTER ... */
acca_first(verbose, out, band, with_shadow,
- count, hist_cold, hist_warm, signa);
+ count, hist_cold, hist_warm, signa);
/* CATEGORIES: NO_DEFINED, WARM_CLOUD, COLD_CLOUD, NULL (= NO_CLOUD) */
value[WARM] = (double)count[WARM] / (double)count[TOTAL];
@@ -101,24 +99,23 @@
value[0] = (double)(count[WARM] + count[COLD]);
idesert = (value[0] == 0. ? 0. : value[0] /
- (value[0] + (double)count[SOIL]));
+ (value[0] + (double)count[SOIL]));
/* BAND-6 CLOUD SIGNATURE DEVELOPMENT */
- if (value[SNOW] > 0.01)
- {
- /* Only the cold clouds are used
- if snow or desert soil is present */
- warm_ambiguous = 1;
+ if (value[SNOW] > 0.01) {
+ /* Only the cold clouds are used
+ if snow or desert soil is present */
+ warm_ambiguous = 1;
}
- else
- {
- /* The cold and warm clouds are combined
- and treated as a single population */
- warm_ambiguous = 0;
- count[COLD] += count[WARM];
- value[COLD] += value[WARM];
- signa[SUM_COLD] += signa[SUM_WARM];
- for (i = 0; i < hist_n; i++) hist_cold[i] += hist_warm[i];
+ else {
+ /* The cold and warm clouds are combined
+ and treated as a single population */
+ warm_ambiguous = 0;
+ count[COLD] += count[WARM];
+ value[COLD] += value[WARM];
+ signa[SUM_COLD] += signa[SUM_WARM];
+ for (i = 0; i < hist_n; i++)
+ hist_cold[i] += hist_warm[i];
}
signa[KMEAN] = (signa[SUM_COLD] / (double)count[COLD]) * SCALE;
@@ -126,81 +123,77 @@
fprintf(stdout, " PRELIMINARY SCENE ANALYSIS\n");
fprintf(stdout, " Desert index: %.3lf\n", idesert);
- fprintf(stdout, " Snow cover : %.3lf %%\n", 100.*value[SNOW]);
- fprintf(stdout, " Cloud cover : %.3lf %%\n", 100.*signa[COVER]);
- fprintf(stdout, " Temperature of %s clouds\n", (warm_ambiguous ? "cold" : "all"));
+ fprintf(stdout, " Snow cover : %.3lf %%\n", 100. * value[SNOW]);
+ fprintf(stdout, " Cloud cover : %.3lf %%\n", 100. * signa[COVER]);
+ fprintf(stdout, " Temperature of %s clouds\n",
+ (warm_ambiguous ? "cold" : "all"));
fprintf(stdout, " Maximum: %.2lf K\n", signa[KMAX]);
fprintf(stdout, " Mean : %.2lf K\n", signa[KMEAN]);
fprintf(stdout, " Minimum: %.2lf K\n", signa[KMIN]);
/* WARNING: variable 'value' reutilizada con nuevos valores */
- if (idesert > 0.5 && signa[COVER] > 0.004 && signa[KMEAN] < 295.)
- {
- value[KUPPER] = quantile( 0.975, hist_cold ) + K_BASE;
- value[KLOWER] = quantile( 0.835, hist_cold ) + K_BASE;
- value[MEAN] = mean(hist_cold) + K_BASE;
- value[VARI] = moment( 2, hist_cold );
- value[SKEW] = moment( 3, hist_cold );
+ if (idesert > 0.5 && signa[COVER] > 0.004 && signa[KMEAN] < 295.) {
+ value[KUPPER] = quantile(0.975, hist_cold) + K_BASE;
+ value[KLOWER] = quantile(0.835, hist_cold) + K_BASE;
+ value[MEAN] = mean(hist_cold) + K_BASE;
+ value[VARI] = moment(2, hist_cold);
+ value[SKEW] = moment(3, hist_cold);
- if (value[SKEW] > 0.)
- {
- shift = value[SKEW];
- if (shift > 1.) shift = 1.;
- shift *= sqrt( value[VARI] );
+ if (value[SKEW] > 0.) {
+ shift = value[SKEW];
+ if (shift > 1.)
+ shift = 1.;
+ shift *= sqrt(value[VARI]);
- x = quantile( 0.9875, hist_cold ) + K_BASE;
- if ((value[KUPPER] + shift) > x)
- {
- value[KUPPER] = x;
- value[KLOWER] = x - shift; /** ??? COMPROBAR ***/
- }
- else
- {
- value[KUPPER] += shift;
- value[KLOWER] += shift;
- }
- }
+ x = quantile(0.9875, hist_cold) + K_BASE;
+ if ((value[KUPPER] + shift) > x) {
+ value[KUPPER] = x;
- fprintf(stdout, " HISTOGRAM CLOUD SIGNATURE\n");
- fprintf(stdout, " Histogram classes: %d\n", hist_n);
- fprintf(stdout, " Mean temperature: %.2lf K\n", value[MEAN]);
- fprintf(stdout, " Standard deviation: %.2lf\n", sqrt(value[VARI]));
- fprintf(stdout, " Skewness: %.2lf\n", value[SKEW]);
- fprintf(stdout, " 97.50 percentile: %.2lf K\n", value[KUPPER]);
- fprintf(stdout, " 83.50 percentile: %.2lf K\n", value[KLOWER]);
-// fprintf(stdout, " 50.00 percentile: %.2lf K\n", quantile(0.5,hist_cold) + K_BASE);
+ value[KLOWER] = x - shift; /** ??? COMPROBAR ***/
+ }
+ else {
+ value[KUPPER] += shift;
+ value[KLOWER] += shift;
+ }
+ }
+
+ fprintf(stdout, " HISTOGRAM CLOUD SIGNATURE\n");
+ fprintf(stdout, " Histogram classes: %d\n", hist_n);
+ fprintf(stdout, " Mean temperature: %.2lf K\n", value[MEAN]);
+ fprintf(stdout, " Standard deviation: %.2lf\n",
+ sqrt(value[VARI]));
+ fprintf(stdout, " Skewness: %.2lf\n", value[SKEW]);
+ fprintf(stdout, " 97.50 percentile: %.2lf K\n", value[KUPPER]);
+ fprintf(stdout, " 83.50 percentile: %.2lf K\n", value[KLOWER]);
+ // fprintf(stdout, " 50.00 percentile: %.2lf K\n", quantile(0.5,hist_cold) + K_BASE);
}
- else
- {
- if (signa[KMEAN] < 295.)
- {
- /* Retained warm and cold clouds */
- fprintf(stdout, " Scene with clouds\n");
- warm_ambiguous = 0;
- value[KUPPER] = 0.;
- value[KLOWER] = 0.;
- }
- else
- {
- /* Retained cold clouds */
- fprintf(stdout, " Scene cloud free\n");
- warm_ambiguous = 1;
- value[KUPPER] = 0.;
- value[KLOWER] = 0.;
- }
+ else {
+ if (signa[KMEAN] < 295.) {
+ /* Retained warm and cold clouds */
+ fprintf(stdout, " Scene with clouds\n");
+ warm_ambiguous = 0;
+ value[KUPPER] = 0.;
+ value[KLOWER] = 0.;
+ }
+ else {
+ /* Retained cold clouds */
+ fprintf(stdout, " Scene cloud free\n");
+ warm_ambiguous = 1;
+ value[KUPPER] = 0.;
+ value[KLOWER] = 0.;
+ }
}
/* SECOND FILTER ... */
/* By-pass two processing by retains warm and cold clouds */
- if (two_pass == 0)
- {
- warm_ambiguous = 0;
- value[KUPPER] = 0.;
- value[KLOWER] = 0.;
+ if (two_pass == 0) {
+ warm_ambiguous = 0;
+ value[KUPPER] = 0.;
+ value[KLOWER] = 0.;
}
acca_second(verbose, out, band[BAND6],
- warm_ambiguous, value[KUPPER], value[KLOWER]);
+ warm_ambiguous, value[KUPPER], value[KLOWER]);
/* CATEGORIES: WARM_CLOUD, COLD_CLOUD, NULL (= NO_CLOUD) */
return;
@@ -208,8 +201,8 @@
void acca_first(int verbose, Gfile * out, Gfile band[],
- int with_shadow,
- int count[], int cold[], int warm[], double stats[])
+ int with_shadow,
+ int count[], int cold[], int warm[], double stats[])
{
int i, row, col, nrows, ncols;
@@ -219,7 +212,7 @@
/* Creation of output file */
out->rast = G_allocate_raster_buf(CELL_TYPE);
if ((out->fd = G_open_raster_new(out->name, CELL_TYPE)) < 0)
- G_fatal_error(_("Could not open <%s>"), out->name);
+ G_fatal_error(_("Could not open <%s>"), out->name);
/* ----- ----- */
fprintf(stdout, "Pass one processing ... \n");
@@ -232,111 +225,112 @@
nrows = G_window_rows();
ncols = G_window_cols();
- for (row = 0; row < nrows; row++)
- {
- if (verbose)
- {
- G_percent(row, nrows, 2);
- }
- for (i = BAND2; i <= BAND6; i++)
- {
- if (G_get_d_raster_row(band[i].fd, band[i].rast, row) < 0)
- G_fatal_error(_("Could not read row from <%s>"), band[i].name);
- }
- for (col = 0; col < ncols; col++)
- {
- code = NO_DEFINED;
- /* Null when null pixel in any band */
- for ( i = BAND2; i <= BAND6; i++ )
- {
- if (G_is_d_null_value((void *)((DCELL *) band[i].rast + col)))
- {
- code = NO_CLOUD;
- break;
- }
- pixel[i] = (double)((DCELL *) band[i].rast)[col];
- }
- /* Determina los pixeles de sombras */
- if (code == NO_DEFINED && with_shadow)
- {
- code = shadow_algorithm(pixel);
- }
- /* Analiza el valor de los pixeles no definidos */
- if (code == NO_DEFINED)
- {
- code = NO_CLOUD;
- count[TOTAL]++;
- nsdi = (pixel[BAND2] - pixel[BAND5]) /
- (pixel[BAND2] + pixel[BAND5]);
- /* ----------------------------------------------------- */
- /* Brightness Threshold: Eliminates dark images */
- if (pixel[BAND3] > th_1)
- {
- /* Normalized Snow Difference Index: Eliminates many types of snow */
- if (nsdi > th_2[0] && nsdi < th_2[1])
- {
- /* Temperature Threshold: Eliminates warm image features */
- if (pixel[BAND6] < th_3)
- {
- rat56 = (1 - pixel[BAND4]) * pixel[BAND6];
- /* Band 5/6 Composite: Eliminates numerous categories including ice */
- if (rat56 < th_4)
- {
- /* Eliminates growing vegetation */
- /* Eliminates senescing vegetation */
- if (pixel[BAND4]/pixel[BAND3] < th_5 &&
- pixel[BAND4]/pixel[BAND2] < th_6)
- {
- rat45 = pixel[BAND4]/pixel[BAND5];
- /* Eliminates rocks and desert */
- if (rat45 > th_7)
- {
- /* Distinguishes warm clouds from cold clouds */
- if (rat56 < th_8)
- {
- code = COLD_CLOUD;
- count[COLD]++;
- /* for statistic */
- stats[SUM_COLD] += (pixel[BAND6]/SCALE);
- hist_put(pixel[BAND6] - K_BASE, cold);
- }
- else
- {
- code = WARM_CLOUD;
- count[WARM]++;
- /* for statistic */
- stats[SUM_WARM] += (pixel[BAND6]/SCALE);
- hist_put(pixel[BAND6] - K_BASE, warm);
- }
- if (pixel[BAND6] > stats[KMAX]) stats[KMAX] = pixel[BAND6];
- if (pixel[BAND6] < stats[KMIN]) stats[KMIN] = pixel[BAND6];
- }
- else { code = NO_DEFINED; count[SOIL]++; }
- }
- else code = NO_DEFINED;
- }
- else code = (pixel[BAND5] < th_4_b) ? NO_CLOUD : NO_DEFINED;
- }
- else code = NO_CLOUD;
- }
- else { code = NO_CLOUD; if (nsdi > th_2_b) count[SNOW]++; }
- }
- else code = (pixel[BAND3] < th_1_b) ? NO_CLOUD : NO_DEFINED;
- /* ----------------------------------------------------- */
- }
- if (code == NO_CLOUD)
- {
- G_set_c_null_value((CELL *) out->rast + col, 1);
- }
- else
- {
- ((CELL *) out->rast)[col] = code;
- }
- }
- if (G_put_raster_row(out->fd, out->rast, CELL_TYPE) < 0)
- {
- G_fatal_error(_("Cannot write row to <%s>"), out->name);
- }
+ for (row = 0; row < nrows; row++) {
+ if (verbose) {
+ G_percent(row, nrows, 2);
+ }
+ for (i = BAND2; i <= BAND6; i++) {
+ if (G_get_d_raster_row(band[i].fd, band[i].rast, row) < 0)
+ G_fatal_error(_("Could not read row from <%s>"),
+ band[i].name);
+ }
+ for (col = 0; col < ncols; col++) {
+ code = NO_DEFINED;
+ /* Null when null pixel in any band */
+ for (i = BAND2; i <= BAND6; i++) {
+ if (G_is_d_null_value((void *)((DCELL *) band[i].rast + col))) {
+ code = NO_CLOUD;
+ break;
+ }
+ pixel[i] = (double)((DCELL *) band[i].rast)[col];
+ }
+ /* Determina los pixeles de sombras */
+ if (code == NO_DEFINED && with_shadow) {
+ code = shadow_algorithm(pixel);
+ }
+ /* Analiza el valor de los pixeles no definidos */
+ if (code == NO_DEFINED) {
+ code = NO_CLOUD;
+ count[TOTAL]++;
+ nsdi = (pixel[BAND2] - pixel[BAND5]) /
+ (pixel[BAND2] + pixel[BAND5]);
+ /* ----------------------------------------------------- */
+ /* Brightness Threshold: Eliminates dark images */
+ if (pixel[BAND3] > th_1) {
+ /* Normalized Snow Difference Index: Eliminates many types of snow */
+ if (nsdi > th_2[0] && nsdi < th_2[1]) {
+ /* Temperature Threshold: Eliminates warm image features */
+ if (pixel[BAND6] < th_3) {
+ rat56 = (1 - pixel[BAND4]) * pixel[BAND6];
+ /* Band 5/6 Composite: Eliminates numerous categories including ice */
+ if (rat56 < th_4) {
+ /* Eliminates growing vegetation */
+ /* Eliminates senescing vegetation */
+ if (pixel[BAND4] / pixel[BAND3] < th_5 &&
+ pixel[BAND4] / pixel[BAND2] < th_6) {
+ rat45 = pixel[BAND4] / pixel[BAND5];
+ /* Eliminates rocks and desert */
+ if (rat45 > th_7) {
+ /* Distinguishes warm clouds from cold clouds */
+ if (rat56 < th_8) {
+ code = COLD_CLOUD;
+ count[COLD]++;
+ /* for statistic */
+ stats[SUM_COLD] +=
+ (pixel[BAND6] / SCALE);
+ hist_put(pixel[BAND6] - K_BASE,
+ cold);
+ }
+ else {
+ code = WARM_CLOUD;
+ count[WARM]++;
+ /* for statistic */
+ stats[SUM_WARM] +=
+ (pixel[BAND6] / SCALE);
+ hist_put(pixel[BAND6] - K_BASE,
+ warm);
+ }
+ if (pixel[BAND6] > stats[KMAX])
+ stats[KMAX] = pixel[BAND6];
+ if (pixel[BAND6] < stats[KMIN])
+ stats[KMIN] = pixel[BAND6];
+ }
+ else {
+ code = NO_DEFINED;
+ count[SOIL]++;
+ }
+ }
+ else
+ code = NO_DEFINED;
+ }
+ else
+ code =
+ (pixel[BAND5] <
+ th_4_b) ? NO_CLOUD : NO_DEFINED;
+ }
+ else
+ code = NO_CLOUD;
+ }
+ else {
+ code = NO_CLOUD;
+ if (nsdi > th_2_b)
+ count[SNOW]++;
+ }
+ }
+ else
+ code = (pixel[BAND3] < th_1_b) ? NO_CLOUD : NO_DEFINED;
+ /* ----------------------------------------------------- */
+ }
+ if (code == NO_CLOUD) {
+ G_set_c_null_value((CELL *) out->rast + col, 1);
+ }
+ else {
+ ((CELL *) out->rast)[col] = code;
+ }
+ }
+ if (G_put_raster_row(out->fd, out->rast, CELL_TYPE) < 0) {
+ G_fatal_error(_("Cannot write row to <%s>"), out->name);
+ }
}
/* ----- ----- */
@@ -348,7 +342,7 @@
void acca_second(int verbose, Gfile * out, Gfile band,
- int warm_ambiguous, double upper, double lower)
+ int warm_ambiguous, double upper, double lower)
{
int row, col, nrows, ncols;
char *mapset;
@@ -360,77 +354,67 @@
/* Open to read */
mapset = G_find_cell2(out->name, "");
if (mapset == NULL)
- G_fatal_error("cell file [%s] not found", out->name);
+ G_fatal_error("cell file [%s] not found", out->name);
out->rast = G_allocate_raster_buf(CELL_TYPE);
if ((out->fd = G_open_cell_old(out->name, mapset)) < 0)
- G_fatal_error("Cannot open cell file [%s]", out->name);
+ G_fatal_error("Cannot open cell file [%s]", out->name);
/* Open to write */
- sprintf(tmp.name, "_%d.BBB", getpid()) ;
+ sprintf(tmp.name, "_%d.BBB", getpid());
tmp.rast = G_allocate_raster_buf(CELL_TYPE);
if ((tmp.fd = G_open_raster_new(tmp.name, CELL_TYPE)) < 0)
- G_fatal_error(_("Could not open <%s>"), tmp.name);
+ G_fatal_error(_("Could not open <%s>"), tmp.name);
if (upper == 0.)
- fprintf(stdout, "Removing ambiguous pixels ... \n");
+ fprintf(stdout, "Removing ambiguous pixels ... \n");
else
- fprintf(stdout, "Pass two processing ... \n");
+ fprintf(stdout, "Pass two processing ... \n");
nrows = G_window_rows();
ncols = G_window_cols();
- for (row = 0; row < nrows; row++)
- {
- if (verbose)
- {
- G_percent(row, nrows, 2);
- }
- if (G_get_d_raster_row(band.fd, band.rast, row) < 0)
- G_fatal_error(_("Could not read from <%s>"), band.name);
- if (G_get_c_raster_row(out->fd, out->rast, row) < 0)
- G_fatal_error(_("Could not read from <%s>"), out->name);
+ for (row = 0; row < nrows; row++) {
+ if (verbose) {
+ G_percent(row, nrows, 2);
+ }
+ if (G_get_d_raster_row(band.fd, band.rast, row) < 0)
+ G_fatal_error(_("Could not read from <%s>"), band.name);
+ if (G_get_c_raster_row(out->fd, out->rast, row) < 0)
+ G_fatal_error(_("Could not read from <%s>"), out->name);
- for (col = 0; col < ncols; col++)
- {
- if (G_is_c_null_value((void *)((CELL *) out->rast + col)))
- {
- G_set_c_null_value((CELL *) tmp.rast + col, 1);
- }
- else
- {
- code = (int)((CELL *) out->rast)[col];
- /* Resolve ambiguous pixels */
- if (code == NO_DEFINED ||
- (code == WARM_CLOUD) && warm_ambiguous == 1)
- {
- temp = (double)((DCELL *) band.rast)[col];
- if (temp > upper)
- {
- G_set_c_null_value((CELL *) tmp.rast + col, 1);
- }
- else
- {
- if (temp < lower)
- ((CELL *) tmp.rast)[col] = IS_COLD_CLOUD;
- else
- ((CELL *) tmp.rast)[col] = IS_WARM_CLOUD;
- }
- }
- else
- /* Join warm (not ambiguous) and cold clouds */
- if (code == COLD_CLOUD ||
- (code == WARM_CLOUD) && warm_ambiguous == 0)
- {
- ((CELL *) tmp.rast)[col] = IS_COLD_CLOUD;
- }
- else
- ((CELL *) tmp.rast)[col] = IS_SHADOW;
- }
- }
- if (G_put_raster_row(tmp.fd, tmp.rast, CELL_TYPE) < 0)
- {
- G_fatal_error(_("Cannot write to <%s>"), tmp.name);
- }
+ for (col = 0; col < ncols; col++) {
+ if (G_is_c_null_value((void *)((CELL *) out->rast + col))) {
+ G_set_c_null_value((CELL *) tmp.rast + col, 1);
+ }
+ else {
+ code = (int)((CELL *) out->rast)[col];
+ /* Resolve ambiguous pixels */
+ if (code == NO_DEFINED ||
+ (code == WARM_CLOUD) && warm_ambiguous == 1) {
+ temp = (double)((DCELL *) band.rast)[col];
+ if (temp > upper) {
+ G_set_c_null_value((CELL *) tmp.rast + col, 1);
+ }
+ else {
+ if (temp < lower)
+ ((CELL *) tmp.rast)[col] = IS_COLD_CLOUD;
+ else
+ ((CELL *) tmp.rast)[col] = IS_WARM_CLOUD;
+ }
+ }
+ else
+ /* Join warm (not ambiguous) and cold clouds */
+ if (code == COLD_CLOUD ||
+ (code == WARM_CLOUD) && warm_ambiguous == 0) {
+ ((CELL *) tmp.rast)[col] = IS_COLD_CLOUD;
+ }
+ else
+ ((CELL *) tmp.rast)[col] = IS_SHADOW;
+ }
+ }
+ if (G_put_raster_row(tmp.fd, tmp.rast, CELL_TYPE) < 0) {
+ G_fatal_error(_("Cannot write to <%s>"), tmp.name);
+ }
}
/* Finalización */
@@ -464,15 +448,12 @@
int shadow_algorithm(double pixel[])
{
-// return NO_DEFINED;
+ // return NO_DEFINED;
- if( pixel[BAND3] < 0.07
- && (1 - pixel[BAND4]) * pixel[BAND6] > 240.
- && pixel[BAND4] / pixel[BAND2] > 1. // Quita agua 1
-// && (pixel[BAND3] - pixel[BAND5]) / (pixel[BAND3] + pixel[BAND5]) < 0.10
- )
- {
- return IS_SHADOW;
+ if (pixel[BAND3] < 0.07 && (1 - pixel[BAND4]) * pixel[BAND6] > 240. && pixel[BAND4] / pixel[BAND2] > 1. // Quita agua 1
+ // && (pixel[BAND3] - pixel[BAND5]) / (pixel[BAND3] + pixel[BAND5]) < 0.10
+ ) {
+ return IS_SHADOW;
}
return NO_DEFINED;
Modified: grass-addons/imagery/i.landsat.acca/local_proto.h
===================================================================
--- grass-addons/imagery/i.landsat.acca/local_proto.h 2010-09-09 07:50:19 UTC (rev 43432)
+++ grass-addons/imagery/i.landsat.acca/local_proto.h 2010-09-09 07:51:43 UTC (rev 43433)
@@ -21,18 +21,18 @@
typedef struct
{
- int fd;
- void * rast;
- char name[1024];
+ int fd;
+ void *rast;
+ char name[1024];
} Gfile;
-void acca_algorithm(int, Gfile *, Gfile [], int, int);
-void acca_first (int, Gfile *, Gfile [], int, int [], int [], int [], double []);
+void acca_algorithm(int, Gfile *, Gfile[], int, int);
+void acca_first(int, Gfile *, Gfile[], int, int[], int[], int[], double[]);
void acca_second(int, Gfile *, Gfile, int, double, double);
-int shadow_algorithm(double []);
+int shadow_algorithm(double[]);
void filter_holes(int, Gfile *);
Modified: grass-addons/imagery/i.landsat.acca/main.c
===================================================================
--- grass-addons/imagery/i.landsat.acca/main.c 2010-09-09 07:50:19 UTC (rev 43432)
+++ grass-addons/imagery/i.landsat.acca/main.c 2010-09-09 07:51:43 UTC (rev 43433)
@@ -46,7 +46,7 @@
* Check a raster name y return fd of open file
*
*----------------------------------------------*/
-int check_raster( char *raster_name )
+int check_raster(char *raster_name)
{
struct Cell_head cellhd;
RASTER_MAP_TYPE map_type;
@@ -55,28 +55,28 @@
mapset = G_find_cell2(raster_name, "");
if (mapset == NULL) {
- G_message("cell file [%s] not found", raster_name);
- return -1;
+ G_message("cell file [%s] not found", raster_name);
+ return -1;
}
if (G_legal_filename(raster_name) < 0) {
- G_message("[%s] is an illegal name", raster_name);
- return -1;
+ G_message("[%s] is an illegal name", raster_name);
+ return -1;
}
if ((raster_fd = G_open_cell_old(raster_name, mapset)) < 0) {
- G_message("Cannot open cell file [%s]", raster_name);
- return -1;
+ G_message("Cannot open cell file [%s]", raster_name);
+ return -1;
}
if (G_get_cellhd(raster_name, mapset, &cellhd) < 0) {
- G_message("Cannot read file header of [%s]", raster_name);
- return -1;
+ G_message("Cannot read file header of [%s]", raster_name);
+ return -1;
}
if (G_set_window(&cellhd) < 0) {
- G_message("Unable to set region");
- return -1;
+ G_message("Unable to set region");
+ return -1;
}
if ((map_type = G_raster_map_type(raster_name, mapset)) != DCELL_TYPE) {
- G_message("Map is not of DCELL_TYPE");
- return -1;
+ G_message("Map is not of DCELL_TYPE");
+ return -1;
}
return raster_fd;
}
@@ -103,14 +103,16 @@
/* initialize module */
module = G_define_module();
- module->description = _("Landsat TM/ETM+ Automatic Cloud Cover Assessment (ACCA)");
+ module->description =
+ _("Landsat TM/ETM+ Automatic Cloud Cover Assessment (ACCA)");
input = G_define_option();
input->key = _("band_prefix");
input->type = TYPE_STRING;
input->required = YES;
input->gisprompt = _("input,cell,raster");
- input->description = _("Base name of the landsat band rasters ([band_prefix].[band_number])");
+ input->description =
+ _("Base name of the landsat band rasters ([band_prefix].[band_number])");
output = G_define_option();
output->key = _("output");
@@ -124,7 +126,8 @@
hist->type = TYPE_INTEGER;
hist->required = NO;
hist->gisprompt = _("input,integer");
- hist->description = _("Number of classes in the cloud temperature histogram");
+ hist->description =
+ _("Number of classes in the cloud temperature histogram");
hist->answer = "100";
sat5 = G_define_flag();
@@ -148,55 +151,52 @@
shadow->answer = 0;
if (G_parser(argc, argv))
- exit(EXIT_FAILURE);
+ exit(EXIT_FAILURE);
/* stores OPTIONS and FLAGS to variables */
- hist_n = atoi( hist->answer );
- if( hist_n < 10) hist_n = 10;
+ hist_n = atoi(hist->answer);
+ if (hist_n < 10)
+ hist_n = 10;
- in_name = input->answer;
+ in_name = input->answer;
- for (i = BAND2; i <= BAND6; i++)
- {
- snprintf(band[i].name, 127, "%s.%d%c", in_name, i+2,
- (i==BAND6 && !sat5->answer ? '1' : '\0'));
- if ((band[i].fd = check_raster(band[i].name)) < 0)
- {
- G_fatal_error(_("Error in filename [%s]!"), band[i].name);
- }
- band[i].rast = G_allocate_raster_buf(DCELL_TYPE);
+ for (i = BAND2; i <= BAND6; i++) {
+ snprintf(band[i].name, 127, "%s.%d%c", in_name, i + 2,
+ (i == BAND6 && !sat5->answer ? '1' : '\0'));
+ if ((band[i].fd = check_raster(band[i].name)) < 0) {
+ G_fatal_error(_("Error in filename [%s]!"), band[i].name);
+ }
+ band[i].rast = G_allocate_raster_buf(DCELL_TYPE);
}
out_name = output->answer;
snprintf(out.name, 127, "%s", out_name);
if (G_legal_filename(out_name) < 0)
- G_fatal_error(_("[%s] is an illegal name"), out.name);
+ G_fatal_error(_("[%s] is an illegal name"), out.name);
/* --------------------------------------- */
-// if( sat5 -> answer )
-// {
-// th_4 = 205.;
-// }
-// acca_test(verbose, &out, band);
- acca_algorithm(verbose, &out, band,
- pass2->answer, shadow->answer);
+ // if( sat5 -> answer )
+ // {
+ // th_4 = 205.;
+ // }
+ // acca_test(verbose, &out, band);
+ acca_algorithm(verbose, &out, band, pass2->answer, shadow->answer);
if (filter->answer)
- filter_holes(verbose, &out);
+ filter_holes(verbose, &out);
/* --------------------------------------- */
- for (i = BAND2; i <= BAND6; i++)
- {
- G_free(band[i].rast);
- G_close_cell(band[i].fd);
+ for (i = BAND2; i <= BAND6; i++) {
+ G_free(band[i].rast);
+ G_close_cell(band[i].fd);
}
-// struct Categories cats;
-// G_read_raster_cats(out.name, char *mapset, cats)
-// G_write_raster_cats(out.name, &cats);
+ // struct Categories cats;
+ // G_read_raster_cats(out.name, char *mapset, cats)
+ // G_write_raster_cats(out.name, &cats);
G_short_history(out.name, "raster", &history);
G_command_history(&history);
Modified: grass-addons/imagery/i.landsat.acca/tools.c
===================================================================
--- grass-addons/imagery/i.landsat.acca/tools.c 2010-09-09 07:50:19 UTC (rev 43432)
+++ grass-addons/imagery/i.landsat.acca/tools.c 2010-09-09 07:51:43 UTC (rev 43433)
@@ -21,17 +21,19 @@
/* Global variable
allow use as parameters of command line */
-int hist_n = 100; /* interval 100/hist_n */
+int hist_n = 100; /* interval 100/hist_n */
void hist_put(double t, int hist[])
{
int i;
/* scale factor */
- i = (int)(t * ((double)hist_n/100.));
+ i = (int)(t * ((double)hist_n / 100.));
- if (i < 1) i = 1;
- if (i > hist_n) i = hist_n;
+ if (i < 1)
+ i = 1;
+ if (i > hist_n)
+ i = hist_n;
hist[i - 1] += 1;
}
@@ -43,16 +45,15 @@
double mean;
total = 0;
- mean = 0.;
- for( i = 0; i < hist_n; i++ )
- {
- total += hist[i];
- mean += (double)(i * hist[i]);
+ mean = 0.;
+ for (i = 0; i < hist_n; i++) {
+ total += hist[i];
+ mean += (double)(i * hist[i]);
}
mean /= (double)total;
/* remove scale factor */
- return (mean / ((double)hist_n/100.));
+ return (mean / ((double)hist_n / 100.));
}
double moment(int n, int hist[])
@@ -61,24 +62,24 @@
double value, hmean, temp;
total = 0;
- hmean = 0.;
- for( i = 0; i < hist_n; i++ )
- {
- total += hist[i];
- hmean += (double)(i * hist[i]);
+ hmean = 0.;
+ for (i = 0; i < hist_n; i++) {
+ total += hist[i];
+ hmean += (double)(i * hist[i]);
}
hmean /= (double)total;
value = 0.;
- for( i = 0; i < hist_n; i++ )
- {
- temp = 1.;
- for( j = 0; j < n; j++ ) temp *= (i - hmean);
- value += temp * (double)hist[i]/(double)total;
+ for (i = 0; i < hist_n; i++) {
+ temp = 1.;
+ for (j = 0; j < n; j++)
+ temp *= (i - hmean);
+ value += temp * (double)hist[i] / (double)total;
}
/* remove scale factor */
- for( j = 0; j < n; j++ ) value /= ((double)hist_n/100.);
+ for (j = 0; j < n; j++)
+ value /= ((double)hist_n / 100.);
return value;
}
@@ -88,25 +89,22 @@
double value, qmax, qmin;
total = 0;
- for( i = 0; i < hist_n; i++ )
- {
- total += hist[i];
+ for (i = 0; i < hist_n; i++) {
+ total += hist[i];
}
qmax = 1.;
- for( i = hist_n - 1; i >= 0; i-- )
- {
- qmin = qmax - (double)hist[i]/(double)total;
- if( q >= qmin )
- {
- value = (q - qmin) / (qmax - qmin) + (i - 1);
- break;
- }
- qmax = qmin;
+ for (i = hist_n - 1; i >= 0; i--) {
+ qmin = qmax - (double)hist[i] / (double)total;
+ if (q >= qmin) {
+ value = (q - qmin) / (qmax - qmin) + (i - 1);
+ break;
+ }
+ qmax = qmin;
}
/* remove scale factor */
- return (value / ((double)hist_n/100.));
+ return (value / ((double)hist_n / 100.));
}
@@ -119,13 +117,14 @@
#define UNDEFINED -1
-int pval(void * rast, int i)
+int pval(void *rast, int i)
{
- void * ptr = (void *)((CELL *)rast + i);
- if( G_is_c_null_value(ptr) )
- return 0;
+ void *ptr = (void *)((CELL *) rast + i);
+
+ if (G_is_c_null_value(ptr))
+ return 0;
else
- return (int)((CELL *) rast)[i];
+ return (int)((CELL *) rast)[i];
}
void filter_holes(int verbose, Gfile * out)
@@ -134,201 +133,187 @@
char *mapset;
void *arast, *brast, *crast;
- int i, pixel[9], cold, warm, shadow, nulo, lim;
+ int i, pixel[9], cold, warm, shadow, nulo, lim;
Gfile tmp;
nrows = G_window_rows();
ncols = G_window_cols();
- if( nrows < 3 || ncols < 3 )
- return;
+ if (nrows < 3 || ncols < 3)
+ return;
/* Open to read */
mapset = G_find_cell2(out->name, "");
if (mapset == NULL)
- G_fatal_error("cell file [%s] not found", out->name);
+ G_fatal_error("cell file [%s] not found", out->name);
arast = G_allocate_raster_buf(CELL_TYPE);
brast = G_allocate_raster_buf(CELL_TYPE);
crast = G_allocate_raster_buf(CELL_TYPE);
if ((out->fd = G_open_cell_old(out->name, mapset)) < 0)
- G_fatal_error("Cannot open cell file [%s]", out->name);
+ G_fatal_error("Cannot open cell file [%s]", out->name);
/* Open to write */
- sprintf(tmp.name, "_%d.BBB", getpid()) ;
+ sprintf(tmp.name, "_%d.BBB", getpid());
tmp.rast = G_allocate_raster_buf(CELL_TYPE);
if ((tmp.fd = G_open_raster_new(tmp.name, CELL_TYPE)) < 0)
- G_fatal_error(_("Could not open <%s>"), tmp.name);
+ G_fatal_error(_("Could not open <%s>"), tmp.name);
fprintf(stdout, "Filling cloud holes ... \n");
- for (row = 0; row < nrows; row++)
- {
- if (verbose)
- {
- G_percent(row, nrows, 2);
- }
- /* Read row values */
- if (row != 0)
- {
- if (G_get_c_raster_row(out->fd, arast, row - 1) < 0)
- G_fatal_error(_("Could not read from <%s>"), out->name);
- }
- if (G_get_c_raster_row(out->fd, brast, row) < 0)
- {
- G_fatal_error(_("Could not read from <%s>"), out->name);
- }
- if (row != (nrows - 1))
- {
- if (G_get_c_raster_row(out->fd, crast, row + 1) < 0)
- G_fatal_error(_("Could not read from <%s>"), out->name);
- }
- /* Analysis of all pixels */
- for (col = 0; col < ncols; col++)
- {
- pixel[0] = pval(brast, col);
- if( pixel[0] == 0 )
- {
- if( row == 0 )
- {
- pixel[1] = -1;
- pixel[2] = -1;
- pixel[3] = -1;
- if( col == 0 )
- {
- pixel[4] = -1;
- pixel[5] = pval(brast, col + 1);
- pixel[6] = -1;
- pixel[7] = pval(crast, col);
- pixel[8] = pval(crast, col + 1);
- }
- else if( col != (ncols - 1) )
- {
- pixel[4] = pval(brast, col - 1);
- pixel[5] = pval(brast, col + 1);
- pixel[6] = pval(crast, col - 1);
- pixel[7] = pval(crast, col);
- pixel[8] = pval(crast, col + 1);
- }
- else
- {
- pixel[4] = pval(brast, col - 1);
- pixel[5] = -1;
- pixel[6] = pval(crast, col - 1);
- pixel[7] = pval(crast, col);
- pixel[8] = -1;
- }
- }
- else if( row != (nrows - 1) )
- {
- if( col == 0 )
- {
- pixel[1] = -1;
- pixel[2] = pval(arast, col);
- pixel[3] = pval(arast, col + 1);
- pixel[4] = -1;
- pixel[5] = pval(brast, col + 1);
- pixel[6] = -1;
- pixel[7] = pval(crast, col);
- pixel[8] = pval(crast, col + 1);
- }
- else if( col != (ncols - 1) )
- {
- pixel[1] = pval(arast, col - 1);
- pixel[2] = pval(arast, col);
- pixel[3] = pval(arast, col + 1);
- pixel[4] = pval(brast, col - 1);
- pixel[5] = pval(brast, col + 1);
- pixel[6] = pval(crast, col - 1);
- pixel[7] = pval(crast, col);
- pixel[8] = pval(crast, col + 1);
- }
- else
- {
- pixel[1] = pval(arast, col - 1);
- pixel[2] = pval(arast, col);
- pixel[3] = -1;
- pixel[4] = pval(brast, col - 1);
- pixel[5] = -1;
- pixel[6] = pval(crast, col - 1);
- pixel[7] = pval(crast, col);
- pixel[8] = -1;
- }
- }
- else
- {
- pixel[6] = -1;
- pixel[7] = -1;
- pixel[8] = -1;
- if( col == 0 )
- {
- pixel[1] = -1;
- pixel[2] = pval(arast, col);
- pixel[3] = pval(arast, col + 1);
- pixel[4] = -1;
- pixel[5] = pval(brast, col + 1);
- }
- else if( col != (ncols - 1) )
- {
- pixel[1] = pval(arast, col - 1);
- pixel[2] = pval(arast, col);
- pixel[3] = pval(arast, col + 1);
- pixel[4] = pval(brast, col - 1);
- pixel[5] = pval(brast, col + 1);
- }
- else
- {
- pixel[1] = pval(arast, col - 1);
- pixel[2] = pval(arast, col);
- pixel[3] = -1;
- pixel[4] = pval(brast, col - 1);
- pixel[5] = -1;
- }
- }
+ for (row = 0; row < nrows; row++) {
+ if (verbose) {
+ G_percent(row, nrows, 2);
+ }
+ /* Read row values */
+ if (row != 0) {
+ if (G_get_c_raster_row(out->fd, arast, row - 1) < 0)
+ G_fatal_error(_("Could not read from <%s>"), out->name);
+ }
+ if (G_get_c_raster_row(out->fd, brast, row) < 0) {
+ G_fatal_error(_("Could not read from <%s>"), out->name);
+ }
+ if (row != (nrows - 1)) {
+ if (G_get_c_raster_row(out->fd, crast, row + 1) < 0)
+ G_fatal_error(_("Could not read from <%s>"), out->name);
+ }
+ /* Analysis of all pixels */
+ for (col = 0; col < ncols; col++) {
+ pixel[0] = pval(brast, col);
+ if (pixel[0] == 0) {
+ if (row == 0) {
+ pixel[1] = -1;
+ pixel[2] = -1;
+ pixel[3] = -1;
+ if (col == 0) {
+ pixel[4] = -1;
+ pixel[5] = pval(brast, col + 1);
+ pixel[6] = -1;
+ pixel[7] = pval(crast, col);
+ pixel[8] = pval(crast, col + 1);
+ }
+ else if (col != (ncols - 1)) {
+ pixel[4] = pval(brast, col - 1);
+ pixel[5] = pval(brast, col + 1);
+ pixel[6] = pval(crast, col - 1);
+ pixel[7] = pval(crast, col);
+ pixel[8] = pval(crast, col + 1);
+ }
+ else {
+ pixel[4] = pval(brast, col - 1);
+ pixel[5] = -1;
+ pixel[6] = pval(crast, col - 1);
+ pixel[7] = pval(crast, col);
+ pixel[8] = -1;
+ }
+ }
+ else if (row != (nrows - 1)) {
+ if (col == 0) {
+ pixel[1] = -1;
+ pixel[2] = pval(arast, col);
+ pixel[3] = pval(arast, col + 1);
+ pixel[4] = -1;
+ pixel[5] = pval(brast, col + 1);
+ pixel[6] = -1;
+ pixel[7] = pval(crast, col);
+ pixel[8] = pval(crast, col + 1);
+ }
+ else if (col != (ncols - 1)) {
+ pixel[1] = pval(arast, col - 1);
+ pixel[2] = pval(arast, col);
+ pixel[3] = pval(arast, col + 1);
+ pixel[4] = pval(brast, col - 1);
+ pixel[5] = pval(brast, col + 1);
+ pixel[6] = pval(crast, col - 1);
+ pixel[7] = pval(crast, col);
+ pixel[8] = pval(crast, col + 1);
+ }
+ else {
+ pixel[1] = pval(arast, col - 1);
+ pixel[2] = pval(arast, col);
+ pixel[3] = -1;
+ pixel[4] = pval(brast, col - 1);
+ pixel[5] = -1;
+ pixel[6] = pval(crast, col - 1);
+ pixel[7] = pval(crast, col);
+ pixel[8] = -1;
+ }
+ }
+ else {
+ pixel[6] = -1;
+ pixel[7] = -1;
+ pixel[8] = -1;
+ if (col == 0) {
+ pixel[1] = -1;
+ pixel[2] = pval(arast, col);
+ pixel[3] = pval(arast, col + 1);
+ pixel[4] = -1;
+ pixel[5] = pval(brast, col + 1);
+ }
+ else if (col != (ncols - 1)) {
+ pixel[1] = pval(arast, col - 1);
+ pixel[2] = pval(arast, col);
+ pixel[3] = pval(arast, col + 1);
+ pixel[4] = pval(brast, col - 1);
+ pixel[5] = pval(brast, col + 1);
+ }
+ else {
+ pixel[1] = pval(arast, col - 1);
+ pixel[2] = pval(arast, col);
+ pixel[3] = -1;
+ pixel[4] = pval(brast, col - 1);
+ pixel[5] = -1;
+ }
+ }
-// pixel[1] = (row == 0 || col == 0) ? -1 : pval(arast, col - 1);
-// pixel[2] = (row == 0) ? -1 : pval(arast, col);
-// pixel[3] = (row == 0 || col == (ncols-1)) ? -1 : pval(arast, col + 1);
-// pixel[4] = (col == 0) ? -1 : pval(brast, col - 1);
-// pixel[5] = (col == (ncols-1)) ? -1 : pval(brast, col + 1);
-// pixel[6] = (row == (nrows-1) || col == 0) ? -1 : pval(crast, col - 1);
-// pixel[7] = (row == (nrows-1)) ? -1 : pval(crast, col);
-// pixel[8] = (row == (nrows-1) || col == (ncols-1)) ? -1 : pval(crast, col + 1);
+ // pixel[1] = (row == 0 || col == 0) ? -1 : pval(arast, col - 1);
+ // pixel[2] = (row == 0) ? -1 : pval(arast, col);
+ // pixel[3] = (row == 0 || col == (ncols-1)) ? -1 : pval(arast, col + 1);
+ // pixel[4] = (col == 0) ? -1 : pval(brast, col - 1);
+ // pixel[5] = (col == (ncols-1)) ? -1 : pval(brast, col + 1);
+ // pixel[6] = (row == (nrows-1) || col == 0) ? -1 : pval(crast, col - 1);
+ // pixel[7] = (row == (nrows-1)) ? -1 : pval(crast, col);
+ // pixel[8] = (row == (nrows-1) || col == (ncols-1)) ? -1 : pval(crast, col + 1);
- cold = warm = shadow = nulo = 0;
- for( i = 1; i < 9; i++ )
- {
- switch( pixel[i] )
- {
- case IS_COLD_CLOUD: cold++; break;
- case IS_WARM_CLOUD: warm++; break;
- case IS_SHADOW : shadow++; break;
- default: nulo++; break;
- }
- }
- lim = (int)(cold + warm + shadow + nulo)/2;
+ cold = warm = shadow = nulo = 0;
+ for (i = 1; i < 9; i++) {
+ switch (pixel[i]) {
+ case IS_COLD_CLOUD:
+ cold++;
+ break;
+ case IS_WARM_CLOUD:
+ warm++;
+ break;
+ case IS_SHADOW:
+ shadow++;
+ break;
+ default:
+ nulo++;
+ break;
+ }
+ }
+ lim = (int)(cold + warm + shadow + nulo) / 2;
- /* Entra pixel[0] = 0 */
- if( nulo < lim )
- {
- if( shadow >= (cold + warm) ) pixel[0] = IS_SHADOW;
- else pixel[0] = (warm > cold) ? IS_WARM_CLOUD : IS_COLD_CLOUD;
- }
- }
- if( pixel[0] != 0 )
- {
- ((CELL *) tmp.rast)[col] = pixel[0];
- }
- else
- {
- G_set_c_null_value((CELL *) tmp.rast + col, 1);
- }
- }
- if (G_put_raster_row(tmp.fd, tmp.rast, CELL_TYPE) < 0)
- {
- G_fatal_error(_("Cannot write to <%s>"), tmp.name);
- }
+ /* Entra pixel[0] = 0 */
+ if (nulo < lim) {
+ if (shadow >= (cold + warm))
+ pixel[0] = IS_SHADOW;
+ else
+ pixel[0] =
+ (warm > cold) ? IS_WARM_CLOUD : IS_COLD_CLOUD;
+ }
+ }
+ if (pixel[0] != 0) {
+ ((CELL *) tmp.rast)[col] = pixel[0];
+ }
+ else {
+ G_set_c_null_value((CELL *) tmp.rast + col, 1);
+ }
+ }
+ if (G_put_raster_row(tmp.fd, tmp.rast, CELL_TYPE) < 0) {
+ G_fatal_error(_("Cannot write to <%s>"), tmp.name);
+ }
}
G_free(arast);
@@ -353,4 +338,3 @@
return;
}
-
More information about the grass-commit
mailing list