[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[gnuastro-commits] master ba81194 3/4: Input kernel no longer mandatory
|
From: |
Mohammad Akhlaghi |
|
Subject: |
[gnuastro-commits] master ba81194 3/4: Input kernel no longer mandatory for NoiseChisel |
|
Date: |
Fri, 2 Mar 2018 19:04:21 -0500 (EST) |
branch: master
commit ba811949085e51d55c26d0b8578bb36350b24daa
Author: Mohammad Akhlaghi <address@hidden>
Commit: Mohammad Akhlaghi <address@hidden>
Input kernel no longer mandatory for NoiseChisel
NoiseChisel no longer needs an input kernel to operate. In such situations,
from now on, a value of `none' can be passed to the `--kernel' option.
During the work, I noticed that the `INFO_NFF' macro defined in
`bin/noisechisel/clumps.c' doesn't do what its name suggests (count
numbers). It actually keeps the sum of pixels above zero, not the number
(as explained in the comments). So the name and comments of this variable
were also corrected.
This issue (to disable convolution), was motivated by a test run and
discussion with Aaron Watkins.
---
NEWS | 9 ++++++---
THANKS | 1 +
bin/noisechisel/clumps.c | 41 +++++++++++++++++++---------------------
bin/noisechisel/main.h | 2 --
bin/noisechisel/noisechisel.c | 42 +++++++++++++++++++++++++----------------
bin/noisechisel/segmentation.c | 3 ++-
bin/noisechisel/sky.c | 43 ++++++++++++++++++++++--------------------
bin/noisechisel/ui.c | 24 ++++++++++++++++-------
bin/noisechisel/ui.h | 7 +++++++
doc/announce-acknowledge.txt | 1 +
doc/gnuastro.texi | 15 +++++++++------
11 files changed, 111 insertions(+), 77 deletions(-)
diff --git a/NEWS b/NEWS
index 0ca8614..d0b5f3a 100644
--- a/NEWS
+++ b/NEWS
@@ -5,13 +5,16 @@ GNU Astronomy Utilities NEWS -*-
outline -*-
** New features
+ MakeCatalog: The new `--mean' and `--median' options will calculate the
+ mean and median pixel value within an object or clump respectively.
+
+ NoiseChisel: a value of `none' to the `--kernel' option will disable
+ convolution.
+
Statistics: the new `--manualbinrange' allows the bins in histograms or
cumulative frequency plots to be set outside the minimum or maximum
values of the dataset.
- MakeCatalog: The new `--mean' and `--median' options will calculate the
- mean and median pixel value within an object or clump respectively.
-
** Removed features
** Changed features
diff --git a/THANKS b/THANKS
index a5ea730..54c09e7 100644
--- a/THANKS
+++ b/THANKS
@@ -52,6 +52,7 @@ support in Gnuastro. The list is ordered alphabetically (by
family name).
Éric Thiébaut address@hidden
Ignacio Trujillo address@hidden
David Valls-Gabaud address@hidden
+ Aaron Watkins address@hidden
Christopher Willmer address@hidden
Sara Yousefi address@hidden
diff --git a/bin/noisechisel/clumps.c b/bin/noisechisel/clumps.c
index 172f112..76d62cc 100644
--- a/bin/noisechisel/clumps.c
+++ b/bin/noisechisel/clumps.c
@@ -693,7 +693,7 @@ enum infocols
{
INFO_X, /* Flux weighted X center col, 0 by C std. */
INFO_Y, /* Flux weighted Y center col. */
- INFO_NFF, /* Number of non-negative pixels (for X,Y).*/
+ INFO_SFF, /* Sum of non-negative pixels (for X,Y). */
INFO_INFLUX, /* Tatal flux within clump. */
INFO_INAREA, /* Tatal area within clump. */
INFO_RIVFLUX, /* Tatal flux within rivers around clump. */
@@ -731,7 +731,7 @@ clumps_get_raw_info(struct clumps_thread_params *cltprm)
info[ lab * INFO_NCOLS + INFO_INFLUX ] += arr[*a];
if( arr[*a]>0.0f )
{
- info[ lab * INFO_NCOLS + INFO_NFF ] += arr[*a];
+ info[ lab * INFO_NCOLS + INFO_SFF ] += arr[*a];
info[ lab * INFO_NCOLS + INFO_X ] += arr[*a] * (*a/dsize[1]);
info[ lab * INFO_NCOLS + INFO_Y ] += arr[*a] * (*a%dsize[1]);
}
@@ -791,17 +791,17 @@ clumps_get_raw_info(struct clumps_thread_params *cltprm)
/* Especially over the undetected regions, it might happen that
none of the pixels were positive. In that case, set the total
area of the clump to zero so it is no longer considered.*/
- if( row[INFO_NFF]==0.0f ) row[INFO_INAREA]=0;
+ if( row[INFO_SFF]==0.0f ) row[INFO_INAREA]=0;
else
{
- coord[0]=GAL_DIMENSION_FLT_TO_INT(row[INFO_X]/row[INFO_NFF]);
- coord[1]=GAL_DIMENSION_FLT_TO_INT(row[INFO_Y]/row[INFO_NFF]);
+ coord[0]=GAL_DIMENSION_FLT_TO_INT(row[INFO_X]/row[INFO_SFF]);
+ coord[1]=GAL_DIMENSION_FLT_TO_INT(row[INFO_Y]/row[INFO_SFF]);
row[INFO_INSTD] = std[ gal_tile_full_id_from_coord(&p->cp.tl,
coord) ];
/* For a check
printf("---------\n");
- printf("\t%f --> %zu\n", row[INFO_Y]/row[INFO_NFF], coord[1]);
- printf("\t%f --> %zu\n", row[INFO_X]/row[INFO_NFF], coord[0]);
+ printf("\t%f --> %zu\n", row[INFO_Y]/row[INFO_SFF], coord[1]);
+ printf("\t%f --> %zu\n", row[INFO_X]/row[INFO_SFF], coord[0]);
printf("%u: (%zu, %zu): %.3f\n", lab, coord[1]+1,
coord[0]+1, row[INFO_INSTD]);
*/
@@ -809,6 +809,7 @@ clumps_get_raw_info(struct clumps_thread_params *cltprm)
}
}
+
/* Clean up. */
free(dinc);
free(ngblabs);
@@ -830,7 +831,6 @@ clumps_make_sn_table(struct clumps_thread_params *cltprm)
int sky0_det1=cltprm->clprm->sky0_det1;
size_t i, ind, counter=0, infodsize[2]={tablen, INFO_NCOLS};
-
/* If there were no initial clumps, then ignore this function. */
if(cltprm->numinitclumps==0) { cltprm->snind=cltprm->sn=NULL; return; }
@@ -1278,8 +1278,8 @@ clumps_true_find_sn_thresh(struct noisechiselparams *p)
/* Set the extension name. */
switch(clprm.step)
{
- case 1: p->clabel->name = "SKY_CLUMPS_ALL"; break;
- case 2: p->clabel->name = "SKY_CLUMPS_FOR_SN"; break;
+ case 1: p->clabel->name = "SKY_CLUMPS_ALL"; break;
+ case 2: p->clabel->name = "SKY_CLUMPS_FOR_SN"; break;
default:
error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s so "
"we can address the issue. The value %d is not valid for "
@@ -1423,20 +1423,17 @@ clumps_det_label_indexs(struct noisechiselparams *p)
int32_t *a, *l, *lf;
gal_data_t *labindexs=gal_data_array_calloc(p->numdetections+1);
- /* Find the area in each detected objects (to see how much space we need
- to allocate). */
+ /* Find the area in each detected object (to see how much space we need
+ to allocate). If blank values are present, an extra check is
+ necessary, so to get faster results when there aren't any blank
+ values, we'll also do a check. */
areas=gal_data_calloc_array(GAL_TYPE_SIZE_T, p->numdetections+1, __func__,
"areas");
- if(gal_blank_present(p->input, 1))
- {
- lf=(l=p->olabel->array)+p->olabel->size; /* Blank pixels have a */
- do if(*l>0) ++areas[*l]; while(++l<lf); /* negative value in int32. */
- }
- else
- {
- lf=(l=p->olabel->array)+p->olabel->size; do ++areas[*l]; while(++l<lf);
- areas[0]=0;
- }
+ lf=(l=p->olabel->array)+p->olabel->size;
+ do
+ if(*l>0) /* Only labeled regions: *l==0 (undetected), *l<0 (blank). */
+ ++areas[*l];
+ while(++l<lf);
/* For a check.
for(i=0;i<p->numdetections+1;++i)
diff --git a/bin/noisechisel/main.h b/bin/noisechisel/main.h
index 003cacc..1f5d759 100644
--- a/bin/noisechisel/main.h
+++ b/bin/noisechisel/main.h
@@ -37,8 +37,6 @@ along with Gnuastro. If not, see
<http://www.gnu.org/licenses/>.
-
-
/* Main program parameters structure */
struct noisechiselparams
{
diff --git a/bin/noisechisel/noisechisel.c b/bin/noisechisel/noisechisel.c
index 2a67232..a9e4f60 100644
--- a/bin/noisechisel/noisechisel.c
+++ b/bin/noisechisel/noisechisel.c
@@ -61,33 +61,43 @@ noisechisel_convolve(struct noisechiselparams *p)
/* Convovle with sharper kernel. */
if(p->conv==NULL)
{
- /* Make the convolved image. */
- if(!p->cp.quiet) gettimeofday(&t1, NULL);
- p->conv = gal_convolve_spatial(tl->tiles, p->kernel, p->cp.numthreads,
- 1, tl->workoverch);
-
- /* Report and write check images if necessary. */
- if(!p->cp.quiet)
+ /* Do the convolution if a kernel was requested. */
+ if(p->kernel)
{
- if(p->widekernel)
- gal_timing_report(&t1, "Convolved with sharper kernel.", 1);
- else
- gal_timing_report(&t1, "Convolved with given kernel.", 1);
+ /* Make the convolved image. */
+ if(!p->cp.quiet) gettimeofday(&t1, NULL);
+ p->conv = gal_convolve_spatial(tl->tiles, p->kernel,
+ p->cp.numthreads, 1, tl->workoverch);
+
+ /* Report and write check images if necessary. */
+ if(!p->cp.quiet)
+ {
+ if(p->widekernel)
+ gal_timing_report(&t1, "Convolved with sharper kernel.", 1);
+ else
+ gal_timing_report(&t1, "Convolved with given kernel.", 1);
+ }
}
+ else
+ p->conv=p->input;
}
/* Set a fixed name for the convolved image (since it will be used in
many check images). */
- if(p->conv->name) free(p->conv->name);
- gal_checkset_allocate_copy( ( p->widekernel
- ? "CONVOLVED-SHARPER"
- : "CONVOLVED" ), &p->conv->name);
+ if(p->conv!=p->input)
+ {
+ if(p->conv->name) free(p->conv->name);
+ gal_checkset_allocate_copy( ( p->widekernel
+ ? "CONVOLVED-SHARPER"
+ : "CONVOLVED" ), &p->conv->name);
+ }
/* Save the convolution step if necessary. */
if(p->detectionname)
{
gal_fits_img_write(p->input, p->detectionname, NULL, PROGRAM_NAME);
- gal_fits_img_write(p->conv, p->detectionname, NULL, PROGRAM_NAME);
+ if(p->input!=p->conv)
+ gal_fits_img_write(p->conv, p->detectionname, NULL, PROGRAM_NAME);
}
/* Convolve with wider kernel (if requested). */
diff --git a/bin/noisechisel/segmentation.c b/bin/noisechisel/segmentation.c
index 3111a5b..c54453b 100644
--- a/bin/noisechisel/segmentation.c
+++ b/bin/noisechisel/segmentation.c
@@ -875,7 +875,8 @@ segmentation(struct noisechiselparams *p)
if(p->segmentationname)
{
gal_fits_img_write(p->input, p->segmentationname, NULL, PROGRAM_NAME);
- gal_fits_img_write(p->conv, p->segmentationname, NULL, PROGRAM_NAME);
+ if(p->input!=p->conv)
+ gal_fits_img_write(p->conv, p->segmentationname, NULL, PROGRAM_NAME);
p->olabel->name="DETECTION_LABELS";
gal_fits_img_write(p->olabel, p->segmentationname, NULL,
PROGRAM_NAME);
diff --git a/bin/noisechisel/sky.c b/bin/noisechisel/sky.c
index b9beee8..b6e3ab5 100644
--- a/bin/noisechisel/sky.c
+++ b/bin/noisechisel/sky.c
@@ -258,25 +258,28 @@ sky_subtract(struct noisechiselparams *p)
/* First subtract the Sky value from the input image. */
GAL_TILE_PARSE_OPERATE(tile, NULL, 0, 0, {*i-=sky[tid];});
- /* Change to the convolved image. */
- tarray=tile->array;
- tblock=tile->block;
- tile->array=gal_tile_block_relative_to_other(tile, p->conv);
- tile->block=p->conv;
-
- /* The threshold is always low. So for the majority of non-NaN
- pixels in the image, the condition above will be true. If we
- come over a NaN pixel, then by definition of NaN, all
- conditionals will fail.
-
- If an image doesn't have any NaN pixels, only the pixels below
- the threshold have to be checked for a NaN which are by
- definition a very small fraction of the total pixels. And if
- there are NaN pixels in the image. */
- GAL_TILE_PARSE_OPERATE(tile, NULL, 0, 0, {*i-=sky[tid];});
-
- /* Revert back to the original block. */
- tile->array=tarray;
- tile->block=tblock;
+ /* Change to the convolved image (if there is any). */
+ if(p->conv!=p->input)
+ {
+ tarray=tile->array;
+ tblock=tile->block;
+ tile->array=gal_tile_block_relative_to_other(tile, p->conv);
+ tile->block=p->conv;
+
+ /* The threshold is always low. So for the majority of non-NaN
+ pixels in the image, the condition above will be true. If we
+ come over a NaN pixel, then by definition of NaN, all
+ conditionals will fail.
+
+ If an image doesn't have any NaN pixels, only the pixels below
+ the threshold have to be checked for a NaN which are by
+ definition a very small fraction of the total pixels. And if
+ there are NaN pixels in the image. */
+ GAL_TILE_PARSE_OPERATE(tile, NULL, 0, 0, {*i-=sky[tid];});
+
+ /* Revert back to the original block. */
+ tile->array=tarray;
+ tile->block=tblock;
+ }
}
}
diff --git a/bin/noisechisel/ui.c b/bin/noisechisel/ui.c
index bb3ee92..e639090 100644
--- a/bin/noisechisel/ui.c
+++ b/bin/noisechisel/ui.c
@@ -250,7 +250,7 @@ ui_read_check_only_options(struct noisechiselparams *p)
" $ info gnuastro \"Input Output options\"");
/* Kernel checks. */
- if(p->kernelname)
+ if(p->kernelname && strcmp(p->kernelname, UI_NO_CONV_KERNEL_NAME))
{
/* Check if it exists. */
gal_checkset_check_file(p->kernelname);
@@ -466,8 +466,13 @@ ui_prepare_kernel(struct noisechiselparams *p)
/* If a kernel file is given, then use it. Otherwise, use the default
kernel. */
if(p->kernelname)
- p->kernel=gal_fits_img_read_kernel(p->kernelname, p->khdu,
- p->cp.minmapsize);
+ {
+ if( strcmp(p->kernelname, UI_NO_CONV_KERNEL_NAME) )
+ p->kernel=gal_fits_img_read_kernel(p->kernelname, p->khdu,
+ p->cp.minmapsize);
+ else
+ p->kernel=NULL;
+ }
else
{
/* Allocate space for the kernel (we don't want to use the statically
@@ -720,9 +725,14 @@ ui_read_check_inputs_setup(int argc, char *argv[],
else
{
if(p->kernelname)
- printf(" - %s: %s (hdu: %s)\n",
- p->widekernelname ? "Sharp Kernel" : "Kernel",
- p->kernelname, p->khdu);
+ {
+ if( strcmp(p->kernelname, UI_NO_CONV_KERNEL_NAME) )
+ printf(" - %s: %s (hdu: %s)\n",
+ p->widekernelname ? "Sharp Kernel" : "Kernel",
+ p->kernelname, p->khdu);
+ else
+ printf(" - No convolution requested.\n");
+ }
else
printf(" - %s: FWHM=2 pixel Gaussian.\n",
p->widekernelname ? "Sharp Kernel" : "Kernel");
@@ -817,7 +827,6 @@ ui_free_report(struct noisechiselparams *p, struct timeval
*t1)
/* Free the allocated datasets. */
gal_data_free(p->sky);
gal_data_free(p->std);
- gal_data_free(p->conv);
gal_data_free(p->wconv);
gal_data_free(p->input);
gal_data_free(p->kernel);
@@ -825,6 +834,7 @@ ui_free_report(struct noisechiselparams *p, struct timeval
*t1)
gal_data_free(p->olabel);
gal_data_free(p->clabel);
gal_data_free(p->widekernel);
+ if(p->conv!=p->input) gal_data_free(p->conv);
/* Clean up the tile structure. */
p->ltl.numchannels=NULL;
diff --git a/bin/noisechisel/ui.h b/bin/noisechisel/ui.h
index d20c2d3..6e33e94 100644
--- a/bin/noisechisel/ui.h
+++ b/bin/noisechisel/ui.h
@@ -30,6 +30,13 @@ along with Gnuastro. If not, see
<http://www.gnu.org/licenses/>.
+/* Macros. */
+#define UI_NO_CONV_KERNEL_NAME "none"
+
+
+
+
+
/* Option groups particular to this program. */
enum program_args_groups
{
diff --git a/doc/announce-acknowledge.txt b/doc/announce-acknowledge.txt
index 5d21881..1cf9d06 100644
--- a/doc/announce-acknowledge.txt
+++ b/doc/announce-acknowledge.txt
@@ -7,4 +7,5 @@ Ole Streicher
Michel Tallon
Juan C. Tello
Éric Thiébaut
+Aaron Watkins
Sara Yousefi
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 03eb629..e1dcca1 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -13745,12 +13745,15 @@ Gnuastro's options with the @option{--help} option.
@item -k STR
@itemx --kernel=STR
File name of kernel to smooth the image before applying the threshold, see
address@hidden kernel}. The first step of NoiseChisel is to
-convolve/smooth the image and use the convolved image in multiple steps
-during the processing. It will be used to define (and later apply) the
-quantile threshold (see @option{--qthresh}). The convolved image is also
-used to define the clumps (see Section 3.2.1 and Figure 8 of
address@hidden://arxiv.org/abs/1505.01664, Akhlaghi and Ichikawa [2015]}).
address@hidden kernel}. If no convolution is needed, give this option a
+value of @option{none}.
+
+The first step of NoiseChisel is to convolve/smooth the image and use the
+convolved image in multiple steps during the processing. It will be used to
+define (and later apply) the quantile threshold (see
address@hidden). The convolved image is also used to define the clumps
+(see Section 3.2.1 and Figure 8 of @url{https://arxiv.org/abs/1505.01664,
+Akhlaghi and Ichikawa [2015]}).
The @option{--kernel} option is not mandatory. If no kernel is provided, a
2D Gaussian profile with a FWHM of 2 pixels truncated at 5 times the FWHM