[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[gnuastro-commits] master 498e1c3: MakeCatalog: correct treating of NaN
From: |
Mohammad Akhlaghi |
Subject: |
[gnuastro-commits] master 498e1c3: MakeCatalog: correct treating of NaN values pixels |
Date: |
Sat, 23 Nov 2019 11:53:02 -0500 (EST) |
branch: master
commit 498e1c38b26722c15e6a0d32c363d7fb4a2b4793
Author: Mohammad Akhlaghi <address@hidden>
Commit: Mohammad Akhlaghi <address@hidden>
MakeCatalog: correct treating of NaN values pixels
Until now, MakeCatalog's `p->hasblank' was mistakenly only set based on the
labels image, not the values image. Also, if there were NaN pixels in the
values image (and not the labels image), then the brightness and magnitude
error columns and the S/N column would become NaN.
With this commit, we now check the Sky value and the Sky standard deviation
value to be NaN before actually using it. It was thus necessary to keep a
count of how many pixels go into each. Also, when estimating `CCOL_SUM_VAR'
and `OCOL_SUM_VAR' (the sum of the variances including the values dataset)
only non-NaN value and STD pixels are used.
This bug was reported by Joseph Putko.
This fixes bug #57293.
---
NEWS | 1 +
bin/mkcatalog/columns.c | 8 ++---
bin/mkcatalog/main.h | 12 ++++---
bin/mkcatalog/parse.c | 85 ++++++++++++++++++++++++++++++++++---------------
bin/mkcatalog/ui.c | 2 +-
5 files changed, 73 insertions(+), 35 deletions(-)
diff --git a/NEWS b/NEWS
index b8c275e..f7f3c65 100644
--- a/NEWS
+++ b/NEWS
@@ -132,6 +132,7 @@ See the end of the file for license conditions.
bug #57164: MakeCatalog crashes when a label isn't in the dataset.
bug #57180: MakeCatalog reporting infinity S/N when --instd isn't an image.
bug #57200: Generated pkgconfig must request wcslib, not wcs.
+ bug #57293: NaN value for brightness-related columns when values have NaN.
diff --git a/bin/mkcatalog/columns.c b/bin/mkcatalog/columns.c
index a2f2443..666f1c3 100644
--- a/bin/mkcatalog/columns.c
+++ b/bin/mkcatalog/columns.c
@@ -2074,12 +2074,12 @@ columns_fill(struct mkcatalog_passparams *pp)
break;
case UI_KEY_SKY:
- ((float *)colarr)[oind] = MKC_RATIO(oi[OCOL_SUMSKY], oi[OCOL_NUM]);
+ ((float *)colarr)[oind] = MKC_RATIO(oi[OCOL_SUMSKY],
oi[OCOL_NUMSKY]);
break;
case UI_KEY_STD:
((float *)colarr)[oind] = sqrt( MKC_RATIO( oi[ OCOL_SUMVAR ],
- oi[ OCOL_NUM ]) );
+ oi[ OCOL_NUMVAR ]) );
break;
case UI_KEY_SEMIMAJOR:
@@ -2323,12 +2323,12 @@ columns_fill(struct mkcatalog_passparams *pp)
case UI_KEY_SKY:
((float *)colarr)[cind] = MKC_RATIO( ci[ CCOL_SUMSKY],
- ci[ CCOL_NUM] );
+ ci[ CCOL_NUMSKY] );
break;
case UI_KEY_STD:
((float *)colarr)[cind] = sqrt( MKC_RATIO( ci[ CCOL_SUMVAR ],
- ci[ CCOL_NUM ] ));
+ ci[ CCOL_NUMVAR ] ));
break;
case UI_KEY_SEMIMAJOR:
diff --git a/bin/mkcatalog/main.h b/bin/mkcatalog/main.h
index 0aa3f72..e281449 100644
--- a/bin/mkcatalog/main.h
+++ b/bin/mkcatalog/main.h
@@ -75,7 +75,7 @@ enum objectcols
OCOL_NUM, /* Number of values used in this object. */
OCOL_NUMXY, /* Number of values in the first two dims. */
OCOL_SUM, /* Sum of (value-sky) in object. */
- OCOL_SUM_VAR, /* Varience of sum (for brightness error). */
+ OCOL_SUM_VAR, /* Variance including values (not just sky). */
OCOL_MEDIAN, /* Median of value in object. */
OCOL_VX, /* Sum of (value-sky) * x. */
OCOL_VY, /* Sum of (value-sky) * y. */
@@ -84,7 +84,9 @@ enum objectcols
OCOL_VYY, /* Sum of (value-sky) * y * y. */
OCOL_VXY, /* Sum of (value-sky) * x * y. */
OCOL_SUMSKY, /* Sum of sky value on this object. */
+ OCOL_NUMSKY, /* Number of sky value on this object. */
OCOL_SUMVAR, /* Sum of sky variance value on this object. */
+ OCOL_NUMVAR, /* Number of sky value on this object. */
OCOL_SUMWHT, /* Sum of positive image pixels. */
OCOL_NUMWHT, /* Number of positive pixels used for wht. */
OCOL_GX, /* Geometric center of object in X. */
@@ -119,7 +121,7 @@ enum clumpcols
CCOL_NUM, /* Number of values used in clump. */
CCOL_NUMXY, /* Number of values only in first two dims. */
CCOL_SUM, /* River subtracted brightness. */
- CCOL_SUM_VAR, /* Variance of sum (for brightness error). */
+ CCOL_SUM_VAR, /* Variance including values (not just sky). */
CCOL_MEDIAN, /* Median of values in clump. */
CCOL_RIV_NUM, /* Num river pixels around this clump. */
CCOL_RIV_SUM, /* Sum of rivers around clump. */
@@ -130,8 +132,10 @@ enum clumpcols
CCOL_VXX, /* Sum of flux*x*x of this clump. */
CCOL_VYY, /* Sum of flux*y*y of this clump. */
CCOL_VXY, /* Sum of flux*x*y of this clump. */
- CCOL_SUMSKY, /* Sum of sky value on this object. */
- CCOL_SUMVAR, /* Sum of sky variance value on this object. */
+ CCOL_SUMSKY, /* Sum of sky value on this clump. */
+ CCOL_NUMSKY, /* Number of sky value on this clump. */
+ CCOL_SUMVAR, /* Sum of sky variance value on this clump. */
+ CCOL_NUMVAR, /* Number of sky variance value on this clump.*/
CCOL_SUMWHT, /* Sum of positive image pixels for wht. */
CCOL_NUMWHT, /* Num of positive image pixels for wht. */
CCOL_GX, /* Geometric center of clump in X. */
diff --git a/bin/mkcatalog/parse.c b/bin/mkcatalog/parse.c
index 1229de9..8ab9831 100644
--- a/bin/mkcatalog/parse.c
+++ b/bin/mkcatalog/parse.c
@@ -441,10 +441,10 @@ parse_objects(struct mkcatalog_passparams *pp)
double *oi=pp->oi;
gal_data_t *xybin=NULL;
size_t *tsize=pp->tile->dsize;
- uint8_t *xybinarr=NULL, *u, *uf;
- float var, sval, *V=NULL, *SK=NULL, *ST=NULL;
+ uint8_t *u, *uf, goodvalue, *xybinarr=NULL;
size_t d, pind=0, increment=0, num_increment=1;
int32_t *O, *OO, *C=NULL, *objarr=p->objects->array;
+ float var, sval, varval, skyval, *V=NULL, *SK=NULL, *ST=NULL;
float *std=p->std?p->std->array:NULL, *sky=p->sky?p->sky->array:NULL;
/* If tile processing isn't necessary, set `tid' to a blank value. */
@@ -563,8 +563,12 @@ parse_objects(struct mkcatalog_passparams *pp)
/* Value related measurements. */
+ goodvalue=0;
if( p->values && !( p->hasblank && isnan(*V) ) )
{
+ /* For the standard-deviation measurements later. */
+ goodvalue=1;
+
/* General flux summations. */
if(xybin) xybinarr[ pind ]=2;
if(oif[ OCOL_NUM ]) oi[ OCOL_NUM ]++;
@@ -608,13 +612,19 @@ parse_objects(struct mkcatalog_passparams *pp)
/* Sky value based measurements. */
- if(p->sky)
- if(oif[ OCOL_SUMSKY ])
- oi[ OCOL_SUMSKY ] += ( pp->st_sky
- ? *SK /* Full array */
- : ( p->sky->size>1
- ? sky[tid] /* Tile */
- : sky[0] ) ); /* Single value */
+ if(p->sky && oif[ OCOL_SUMSKY ])
+ {
+ skyval = ( pp->st_sky
+ ? (isnan(*SK)?0:*SK) /* Full array
*/
+ : ( p->sky->size>1
+ ? (isnan(sky[tid])?0:sky[tid]) /* Tile
*/
+ : sky[0] ) ); /* Single
value*/
+ if(!isnan(skyval))
+ {
+ oi[ OCOL_NUMSKY ]++;
+ oi[ OCOL_SUMSKY ] += skyval;
+ }
+ }
/* Sky standard deviation based measurements.*/
@@ -622,7 +632,11 @@ parse_objects(struct mkcatalog_passparams *pp)
{
sval = pp->st_std ? *ST : (p->std->size>1?std[tid]:std[0]);
var = p->variance ? sval : sval*sval;
- if(oif[ OCOL_SUMVAR ]) oi[ OCOL_SUMVAR ] += var;
+ if(oif[ OCOL_SUMVAR ] && (!isnan(var)))
+ {
+ oi[ OCOL_NUMVAR ]++;
+ oi[ OCOL_SUMVAR ] += var;
+ }
/* For each pixel, we have a sky contribution to the
counts and the signal's contribution. The standard
deviation in the sky is simply `sval', but the
@@ -633,9 +647,11 @@ parse_objects(struct mkcatalog_passparams *pp)
absolute value, because especially as the signal gets
noisy there will be negative values, and we don't want
them to decrease the variance. */
- if(oif[ OCOL_SUM_VAR ])
- oi[ OCOL_SUM_VAR ] += ( (p->variance ? var : sval)
- + fabs(*V) );
+ if(oif[ OCOL_SUM_VAR ] && goodvalue)
+ {
+ varval=p->variance ? var : sval;
+ if(!isnan(varval)) oi[ OCOL_SUM_VAR ] += varval +
fabs(*V);
+ }
}
}
@@ -718,10 +734,10 @@ parse_clumps(struct mkcatalog_passparams *pp)
gal_data_t *xybin=NULL;
size_t *tsize=pp->tile->dsize;
int32_t *O, *OO, *C=NULL, nlab;
- uint8_t *u, *uf, *cif=p->ciflag;
- float var, sval, *V=NULL, *SK=NULL, *ST=NULL;
+ uint8_t *u, *uf, goodvalue, *cif=p->ciflag;
size_t nngb=gal_dimension_num_neighbors(ndim);
size_t i, ii, d, pind=0, increment=0, num_increment=1;
+ float var, sval, varval, skyval, *V=NULL, *SK=NULL, *ST=NULL;
int32_t *objects=p->objects->array, *clumps=p->clumps->array;
float *std=p->std?p->std->array:NULL, *sky=p->sky?p->sky->array:NULL;
@@ -855,8 +871,12 @@ parse_clumps(struct mkcatalog_passparams *pp)
/* Value related measurements, see `parse_objects' for
comments. */
+ goodvalue=0;
if( p->values && !( p->hasblank && isnan(*V) ) )
{
+ /* For the standard-deviation measurement. */
+ goodvalue=1;
+
/* Fill in the necessary information. */
if(cif[ CCOL_NUM ]) ci[ CCOL_NUM ]++;
if(cif[ CCOL_SUM ]) ci[ CCOL_SUM ] += *V;
@@ -882,13 +902,19 @@ parse_clumps(struct mkcatalog_passparams *pp)
}
/* Sky based measurements. */
- if(p->sky)
- if(cif[ CCOL_SUMSKY ])
- ci[ CCOL_SUMSKY ] += ( pp->st_sky
- ? *SK /* Full */
- : ( p->sky->size>1
- ? sky[tid] /* Tile */
- : sky[0] ) ); /* 1 value */
+ if(p->sky && cif[ CCOL_SUMSKY ])
+ {
+ skyval = ( pp->st_sky
+ ? *SK /* Full. */
+ : ( p->sky->size>1
+ ? sky[tid] /* Tile. */
+ : sky[0] ) ); /* 1 value. */
+ if(!isnan(skyval))
+ {
+ ci[ CCOL_NUMSKY ]++;
+ ci[ CCOL_SUMSKY ] += skyval;
+ }
+ }
/* Sky Standard deviation based measurements, see
`parse_objects' for comments. */
@@ -898,10 +924,17 @@ parse_clumps(struct mkcatalog_passparams *pp)
? *ST
: (p->std->size>1 ? std[tid] : std[0]) );
var = p->variance ? sval : sval*sval;
- if(cif[ CCOL_SUMVAR ]) ci[ CCOL_SUMVAR ] += var;
- if(cif[ CCOL_SUM_VAR ])
- ci[ CCOL_SUM_VAR ] += ( (p->variance ? var : sval)
- + fabs(*V) );
+ if(cif[ CCOL_SUMVAR ] && (!isnan(var)))
+ {
+ ci[ CCOL_NUMVAR ]++;
+ ci[ CCOL_SUMVAR ] += var;
+ }
+ if(cif[ CCOL_SUM_VAR ] && goodvalue)
+ {
+ varval=p->variance ? var : sval;
+ if(!isnan(varval))
+ ci[ CCOL_SUM_VAR ] += varval + fabs(*V);
+ }
}
}
diff --git a/bin/mkcatalog/ui.c b/bin/mkcatalog/ui.c
index adcfe46..c6c1a55 100644
--- a/bin/mkcatalog/ui.c
+++ b/bin/mkcatalog/ui.c
@@ -1159,7 +1159,7 @@ ui_preparations_read_inputs(struct mkcatalogparams *p)
/* Initially, `p->hasblank' was set based on the objects image, but
it may happen that the objects image only has zero values for
blank pixels, so we'll also do a check on the input image. */
- p->hasblank = gal_blank_present(p->objects, 1);
+ p->hasblank = gal_blank_present(p->values, 1);
/* Reset the units of the value-based columns if the input dataset
has defined units. */
[Prev in Thread] |
Current Thread |
[Next in Thread] |
- [gnuastro-commits] master 498e1c3: MakeCatalog: correct treating of NaN values pixels,
Mohammad Akhlaghi <=