[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[gnuastro-commits] master e61e3526 1/3: Library (convolve.h): new argume
|
From: |
Mohammad Akhlaghi |
|
Subject: |
[gnuastro-commits] master e61e3526 1/3: Library (convolve.h): new argument for convolution over blanks |
|
Date: |
Wed, 18 Oct 2023 09:53:54 -0400 (EDT) |
branch: master
commit e61e352615761fd17ee5df1b46a119241d31bf50
Author: Mohammad Akhlaghi <mohammad@akhlaghi.org>
Commit: Mohammad Akhlaghi <mohammad@akhlaghi.org>
Library (convolve.h): new argument for convolution over blanks
Until now, the library's convolution functions would always ignore blank
pixels. In most scenarios this is fine and expected, but in some
situations, users may want to use convolution as a proxy for interpolation
or dilation.
With this commit, a new 'conv_on_blank' argument has been added to the core
convolution library to allow for this feature when the user wants this
behavior. Also, the Convolve program now has a new '--conv-on-blank' option
for users to directly activate this feature on the command-line.
This was suggested by Raul Infante-Sainz.
---
NEWS | 10 ++++++++++
bin/convolve/args.h | 13 +++++++++++++
bin/convolve/convolve.c | 3 ++-
bin/convolve/main.h | 1 +
bin/convolve/ui.h | 1 +
bin/noisechisel/noisechisel.c | 6 ++++--
bin/segment/segment.c | 3 ++-
bin/statistics/sky.c | 3 ++-
doc/gnuastro.texi | 16 +++++++++++++---
lib/convolve.c | 22 ++++++++++++++--------
lib/gnuastro/convolve.h | 4 +++-
lib/tile.c | 4 ++--
12 files changed, 67 insertions(+), 19 deletions(-)
diff --git a/NEWS b/NEWS
index 5eef5683..01312391 100644
--- a/NEWS
+++ b/NEWS
@@ -96,6 +96,12 @@ See the end of the file for license conditions.
Gnuastro could only read TIFF files). This step was written by Fathma
Mehnoor.
+*** Convolve
+ --conv-on-blank: do not ignore blank pixels in the convolution. This will
+ effectively expand the non-blank regions of your dataset into the blank
+ regions and only works in spatial-domain convolution. This was
+ suggested by Raul Infante-Sainz.
+
*** MakeCatalog
- New measurements:
--river-min: minimum river value around a clump.
@@ -229,6 +235,10 @@ See the end of the file for license conditions.
Implemented by Sepideh Eskandarlou
*** Library
+ - gal_convolve_spatial: new argument to allow convolution over blank
+ elements.
+ - gal_convolve_spatial_correct_ch_edge: new argument to allow convolution
+ over blank elements.
- gal_polygon_area_flat: new name for the old 'gal_polygon_area'
function. This was necessary because this function assumes flat
coordinates and can lead to wrong results when used with a polygon
diff --git a/bin/convolve/args.h b/bin/convolve/args.h
index 04c4d605..ce23c6bf 100644
--- a/bin/convolve/args.h
+++ b/bin/convolve/args.h
@@ -155,6 +155,19 @@ struct argp_option program_options[] =
GAL_OPTIONS_NOT_MANDATORY,
GAL_OPTIONS_NOT_SET
},
+ {
+ "conv-on-blank",
+ UI_KEY_CONVONBLANK,
+ 0,
+ 0,
+ "Do convolution on blank pixels also.",
+ GAL_OPTIONS_GROUP_OUTPUT,
+ &p->conv_on_blank,
+ GAL_OPTIONS_NO_ARG_TYPE,
+ GAL_OPTIONS_RANGE_0_OR_1,
+ GAL_OPTIONS_NOT_MANDATORY,
+ GAL_OPTIONS_NOT_SET
+ },
diff --git a/bin/convolve/convolve.c b/bin/convolve/convolve.c
index ad98d511..8cd44ef4 100644
--- a/bin/convolve/convolve.c
+++ b/bin/convolve/convolve.c
@@ -800,7 +800,8 @@ convolve(struct convolveparams *p)
out=gal_convolve_spatial(multidim ? cp->tl.tiles : p->input, p->kernel,
cp->numthreads,
multidim ? !p->noedgecorrection : 1,
- multidim ? cp->tl.workoverch : 1 );
+ multidim ? cp->tl.workoverch : 1,
+ p->conv_on_blank);
/* Clean up: free the actual input and replace it's pointer with the
convolved dataset to save as output. */
diff --git a/bin/convolve/main.h b/bin/convolve/main.h
index bbab8221..370f50e3 100644
--- a/bin/convolve/main.h
+++ b/bin/convolve/main.h
@@ -83,6 +83,7 @@ struct convolveparams
char *domainstr; /* String value specifying domain. */
size_t makekernel; /* Make a kernel to create input. */
uint8_t noedgecorrection; /* Do not correct spatial edge effects. */
+ uint8_t conv_on_blank; /* Do convolution on blank pixels also. */
/* Internal */
int isfits; /* Input is a FITS file. */
diff --git a/bin/convolve/ui.h b/bin/convolve/ui.h
index 888f93c8..b70c7b05 100644
--- a/bin/convolve/ui.h
+++ b/bin/convolve/ui.h
@@ -62,6 +62,7 @@ enum option_keys_enum
/* Only with long version (start with a value 1000, the rest will be set
automatically). */
UI_KEY_KERNELCOLUMN = 1000,
+ UI_KEY_CONVONBLANK,
UI_KEY_NOKERNELFLIP,
UI_KEY_NOKERNELNORM,
UI_KEY_NOEDGECORRECTION,
diff --git a/bin/noisechisel/noisechisel.c b/bin/noisechisel/noisechisel.c
index fa9456e5..5fe4d759 100644
--- a/bin/noisechisel/noisechisel.c
+++ b/bin/noisechisel/noisechisel.c
@@ -68,7 +68,8 @@ noisechisel_convolve(struct noisechiselparams *p)
/* 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);
+ p->cp.numthreads, 1,
+ tl->workoverch, 0);
/* Report and write check images if necessary. */
if(!p->cp.quiet)
@@ -106,7 +107,8 @@ noisechisel_convolve(struct noisechiselparams *p)
{
if(!p->cp.quiet) gettimeofday(&t1, NULL);
p->wconv=gal_convolve_spatial(tl->tiles, p->widekernel,
- p->cp.numthreads, 1, tl->workoverch);
+ p->cp.numthreads, 1,
+ tl->workoverch, 0);
gal_checkset_allocate_copy("CONVOLVED-WIDER", &p->wconv->name);
if(!p->cp.quiet)
diff --git a/bin/segment/segment.c b/bin/segment/segment.c
index fd49eab0..6b32f646 100644
--- a/bin/segment/segment.c
+++ b/bin/segment/segment.c
@@ -75,7 +75,8 @@ segment_convolve(struct segmentparams *p)
/* 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);
+ p->cp.numthreads, 1,
+ tl->workoverch, 0);
/* Report and write check images if necessary. */
if(!p->cp.quiet)
diff --git a/bin/statistics/sky.c b/bin/statistics/sky.c
index 992c6ac2..cc2ac1f0 100644
--- a/bin/statistics/sky.c
+++ b/bin/statistics/sky.c
@@ -170,7 +170,8 @@ sky(struct statisticsparams *p)
{
if(!cp->quiet) gettimeofday(&t1, NULL);
p->convolved=gal_convolve_spatial(tl->tiles, p->kernel,
- cp->numthreads, 1, tl->workoverch);
+ cp->numthreads, 1,
+ tl->workoverch, 0);
if(p->checksky)
gal_fits_img_write(p->convolved, p->checkskyname, NULL,
PROGRAM_NAME);
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index e9e4600f..9f2e08dd 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -23686,6 +23686,15 @@ By default, the input kernel is flipped to avoid the
output getting flipped; see
Do not normalize the kernel after reading it, such that the sum of its pixels
is unity.
As described in @ref{Convolution process}, the kernel is normalized by default.
+@item --conv-on-blank
+Do not ignore blank pixels in the convolution.
+The output pixels that were originally non-blank are not affected by this
option (they will have the same value if this option is called or not).
+This option just expands/dilates the non-blank regions of your dataset into
the blank regions and only works in spatial domain convolution.
+Therefore, with this option convolution can be used as a proxy for
interpolation or dilation.
+
+By default, blank pixels are ignored during spatial domain convolution; so the
input and output have exactly the same number of blank pixels.
+With this option, the blank pixels that are sufficiently close to non-blank
pixels (based on the kernel) will be given a value based on the non-blank
elements that overlap with the kernel for that blank pixel (see @ref{Edges in
the spatial domain}).
+
@item -d STR
@itemx --domain=STR
@cindex Discrete Fourier transform
@@ -23695,7 +23704,6 @@ The acceptable values are `@code{spatial}' and
`@code{frequency}', corresponding
For large images, the frequency domain process will be more efficient than
convolving in the spatial domain.
However, the edges of the image will loose some flux (see @ref{Edges in the
spatial domain}) and the image must not contain any blank pixels, see
@ref{Spatial vs. Frequency domain}.
-
@item --checkfreqsteps
With this option a file with the initial name of the output file will be
created that is suffixed with @file{_freqsteps.fits}, all the steps done to
arrive at the final convolved image are saved as extensions in this file.
The extensions in order are:
@@ -41408,9 +41416,10 @@ Because of the complete introduction that was
presented there, we will directly
As of this version, only spatial domain convolution is available in Gnuastro's
libraries.
We have not had the time to liberate the frequency domain function convolution
and de-convolution functions that are available in the Convolve
program@footnote{Hence any help would be greatly appreciated.}.
-@deftypefun {gal_data_t *} gal_convolve_spatial (gal_data_t @code{*tiles},
gal_data_t @code{*kernel}, size_t @code{numthreads}, int @code{edgecorrection},
int @code{convoverch})
+@deftypefun {gal_data_t *} gal_convolve_spatial (gal_data_t @code{*tiles},
gal_data_t @code{*kernel}, size_t @code{numthreads}, int @code{edgecorrection},
int @code{convoverch}, int @code{conv_on_blank})
Convolve the given @code{tiles} dataset (possibly a list of tiles, see
@ref{List of gal_data_t} and @ref{Tessellation library}) with @code{kernel} on
@code{numthreads} threads.
When @code{edgecorrection} is non-zero, it will correct for the edge dimming
effects as discussed in @ref{Edges in the spatial domain}.
+When @code{conv_on_blank} is non-zero, this function will also attempt
convolution over the blank pixels (and therefore give values to the blank
pixels that are near non-blank pixels).
@code{tiles} can be a single/complete dataset, but in that case the speed will
be very slow.
Therefore, for larger images, it is recommended to give a list of tiles
covering a dataset.
@@ -41425,13 +41434,14 @@ This behavior may be disabled when @code{convoverch}
is non-zero.
In this case, it will ignore channel borders (if they exist) and mix all
pixels that cover the kernel within the dataset.
@end deftypefun
-@deftypefun void gal_convolve_spatial_correct_ch_edge (gal_data_t
@code{*tiles}, gal_data_t @code{*kernel}, size_t @code{numthreads}, int
@code{edgecorrection}, gal_data_t @code{*tocorrect})
+@deftypefun void gal_convolve_spatial_correct_ch_edge (gal_data_t
@code{*tiles}, gal_data_t @code{*kernel}, size_t @code{numthreads}, int
@code{edgecorrection}, int @code{conv_on_blank}, gal_data_t @code{*tocorrect})
Correct the edges of channels in an already convolved image when it was
initially convolved with @code{gal_convolve_spatial} and @code{convoverch==0}.
In that case, strong boundaries might exist on the channel edges.
So if you later need to remove those boundaries at later steps of your
processing, you can call this function.
It will only do convolution on the tiles that are near the edge and were
effected by the channel borders.
Other pixels in the image will not be touched.
Hence, it is much faster.
+When @code{conv_on_blank} is non-zero, this function will also attempt
convolution over the blank pixels (and therefore give values to the blank
pixels that are near non-blank pixels).
@end deftypefun
@node Pooling functions, Interpolation, Convolution functions, Gnuastro library
diff --git a/lib/convolve.c b/lib/convolve.c
index 02074f8d..8e44d662 100644
--- a/lib/convolve.c
+++ b/lib/convolve.c
@@ -99,8 +99,8 @@ struct per_thread_spatial_prm
size_t *overlap_start; /* Starting coordinate of kernel overlap. */
size_t *kernel_start; /* Kernel starting point. */
size_t *host_start; /* Starting coordinate of host. */
- size_t *pix; /* 2*ndim: starting and ending of tile,
- Later, just the pixel being convolved. */
+ size_t *pix; /* 2*ndim: starting and ending of tile, */
+ /* Later, just the pixel being convolved. */
int on_edge; /* If the tile is on the edge or not. */
gal_data_t *host; /* Size of host (channel or block). */
struct spatial_params *cprm; /* Link to main structure for all threads. */
@@ -120,6 +120,7 @@ struct spatial_params
gal_data_t *tocorrect; /* (possible) convolved image to correct. */
int convoverch; /* Ignore channel edges in convolution. */
int edgecorrection; /* Correct convolution's edge effects. */
+ uint8_t conv_on_blank; /* Do convolution over blank pixels also. */
struct per_thread_spatial_prm *pprm; /* Array of per-thread parameters. */
};
@@ -358,7 +359,7 @@ convolve_spatial_tile(struct per_thread_spatial_prm *pprm)
/* If the input on this pixel is a NaN, then just set the output
to NaN too and go onto the next pixel. 'in_v' is the pointer
on this pixel. */
- if( isnan(*in_v) )
+ if( isnan(*in_v) && cprm->conv_on_blank==0 )
out[ in_v - in ]=NAN;
else
{
@@ -487,7 +488,8 @@ convolve_spatial_on_thread(void *inparam)
static gal_data_t *
gal_convolve_spatial_general(gal_data_t *tiles, gal_data_t *kernel,
size_t numthreads, int edgecorrection,
- int convoverch, gal_data_t *tocorrect)
+ int convoverch, uint8_t conv_on_blank,
+ gal_data_t *tocorrect)
{
struct spatial_params params;
gal_data_t *out, *block=gal_tile_block(tiles);
@@ -536,6 +538,7 @@ gal_convolve_spatial_general(gal_data_t *tiles, gal_data_t
*kernel,
params.kernel=kernel;
params.tocorrect=tocorrect;
params.convoverch=convoverch;
+ params.conv_on_blank=conv_on_blank;
params.edgecorrection=edgecorrection;
@@ -569,7 +572,8 @@ gal_convolve_spatial_general(gal_data_t *tiles, gal_data_t
*kernel,
input, the 'next' element has to be 'NULL'. */
gal_data_t *
gal_convolve_spatial(gal_data_t *tiles, gal_data_t *kernel,
- size_t numthreads, int edgecorrection, int convoverch)
+ size_t numthreads, int edgecorrection, int convoverch,
+ int conv_on_blank)
{
/* When there isn't any tile structure, 'convoverch' must be set to
one. Recall that the input can be a single full dataset also. */
@@ -577,7 +581,8 @@ gal_convolve_spatial(gal_data_t *tiles, gal_data_t *kernel,
/* Call the general function. */
return gal_convolve_spatial_general(tiles, kernel, numthreads,
- edgecorrection, convoverch, NULL);
+ edgecorrection, convoverch,
+ conv_on_blank, NULL);
}
@@ -593,7 +598,7 @@ gal_convolve_spatial(gal_data_t *tiles, gal_data_t *kernel,
void
gal_convolve_spatial_correct_ch_edge(gal_data_t *tiles, gal_data_t *kernel,
size_t numthreads, int edgecorrection,
- gal_data_t *tocorrect)
+ int conv_on_blank, gal_data_t *tocorrect)
{
gal_data_t *block=gal_tile_block(tiles);
@@ -609,5 +614,6 @@ gal_convolve_spatial_correct_ch_edge(gal_data_t *tiles,
gal_data_t *kernel,
/* Call the general function, which will do the correction. */
gal_convolve_spatial_general(tiles, kernel, numthreads,
- edgecorrection, 0, tocorrect);
+ edgecorrection, 0, conv_on_blank,
+ tocorrect);
}
diff --git a/lib/gnuastro/convolve.h b/lib/gnuastro/convolve.h
index 9d298670..aa64d759 100644
--- a/lib/gnuastro/convolve.h
+++ b/lib/gnuastro/convolve.h
@@ -46,12 +46,14 @@ __BEGIN_C_DECLS /* From C++ preparations */
gal_data_t *
gal_convolve_spatial(gal_data_t *tiles, gal_data_t *kernel,
- size_t numthreads, int edgecorrection, int convoverch);
+ size_t numthreads, int edgecorrection,
+ int convoverch, int conv_on_blank);
void
gal_convolve_spatial_correct_ch_edge(gal_data_t *tiles, gal_data_t *kernel,
size_t numthreads, int edgecorrection,
+ int conv_on_blank,
gal_data_t *tocorrect);
diff --git a/lib/tile.c b/lib/tile.c
index 7dc67b48..6c5eda2b 100644
--- a/lib/tile.c
+++ b/lib/tile.c
@@ -1185,7 +1185,7 @@ gal_tile_full_values_smooth(gal_data_t *tilevalues,
/* Do the smoothing. */
if(tl->workoverch)
- smoothed=gal_convolve_spatial(tilevalues, kernel, numthreads, 1, 1);
+ smoothed=gal_convolve_spatial(tilevalues, kernel, numthreads, 1, 1, 0);
else
{
/* Create the tile structure. */
@@ -1195,7 +1195,7 @@ gal_tile_full_values_smooth(gal_data_t *tilevalues,
gal_tile_full_two_layers(tilevalues, &ttl);
/* Do the convolution separately on each channel. */
- smoothed=gal_convolve_spatial(ttl.tiles, kernel, numthreads, 1, 0);
+ smoothed=gal_convolve_spatial(ttl.tiles, kernel, numthreads, 1, 0, 0);
/* Clean up. */
ttl.tilesize=ttl.numchannels=NULL;