[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[gnuastro-commits] master b4c2b5a: gal_box_bound_ellipse return width is
From: |
Mohammad Akhlaghi |
Subject: |
[gnuastro-commits] master b4c2b5a: gal_box_bound_ellipse return width isn't extra |
Date: |
Mon, 23 Apr 2018 09:50:46 -0400 (EDT) |
branch: master
commit b4c2b5a9099564c844b5ad74257675266cdab2c1
Author: Mohammad Akhlaghi <address@hidden>
Commit: Mohammad Akhlaghi <address@hidden>
gal_box_bound_ellipse return width isn't extra
Until now, the width returned by `gal_box_bound_ellipse' would be 2 pixels
extra along each dimension, thus leaving some extra pixels on the edges
when MakeProfiles used it. This was done to accommodate extreme cases,
where the center of a profile was just on the edge of the pixel and the
trunctation radius also had a fractional value.
However, after some checks, we see that even in the most extreme cases,
only a single (very distant) elliptical layer may be lost (not at all
significant). But in most (normal cases), these extra layers will cause
many extra pixels which can be a burden on MakeProfiles its self or
higher-level programs (for example when building a kernel to convolve
with).
Also, while checking the `gal_box_bound_ellipse' functions, I noticed that
in the documentation, we were mistakenly describing `a' and `b' as
major/minor axis length. While in practice they are `semi-major' and
`semi-minor'.
Since there is no extra layer, we also don't need the call to Crop when
making NoiseChisel's default internal kernel (in
`lib/gnuastro-internal/kernel-2d.h'). A new kernel is now created and used
(with a larger number of random points).
---
bin/mkprof/oneprofile.c | 23 +++++++++++++++++-
doc/gnuastro.texi | 27 +++++++++++----------
lib/box.c | 4 ++--
lib/gnuastro-internal/kernel-2d.h | 50 +++++++++++++++++++--------------------
4 files changed, 63 insertions(+), 41 deletions(-)
diff --git a/bin/mkprof/oneprofile.c b/bin/mkprof/oneprofile.c
index 865e63f..aab2524 100644
--- a/bin/mkprof/oneprofile.c
+++ b/bin/mkprof/oneprofile.c
@@ -172,7 +172,9 @@ oneprofile_r_circle(size_t index, struct mkonthread *mkp)
float
oneprofile_randompoints(struct mkonthread *mkp)
{
+ double r_before=mkp->r;
double range[2], sum=0.0f;
+ double coord_before[2]={mkp->coord[0], mkp->coord[1]};
size_t i, j, numrandom=mkp->p->numrandom, ndim=mkp->p->ndim;
/* Set the range in each dimension. */
@@ -188,6 +190,15 @@ oneprofile_randompoints(struct mkonthread *mkp)
sum+=mkp->profile(mkp);
}
+ /* Reset the original distance and coordinate of the pixel and return the
+ average random value. The resetting is mostly redundant (only useful
+ in checks), but since it has a very negligible cost (compared to the
+ random checks above) cost, its good to reset it to help in debugging
+ when necessary (avoid confusion when un-commenting the checks in
+ `oneprofile_pix_by_pix'). */
+ mkp->r=r_before;
+ mkp->coord[0]=coord_before[0];
+ mkp->coord[1]=coord_before[1];
return sum/numrandom;
}
@@ -376,7 +387,7 @@ oneprofile_pix_by_pix(struct mkonthread *mkp)
oneprofile_r_el(mkp);
if(mkp->r > truncr) continue;
- /* Find the value for this pixel: */
+ /* Set the range for this pixel. */
for(i=0;i<ndim;++i)
{
mkp->lower[i] = mkp->coord[i] - hp;
@@ -389,6 +400,12 @@ oneprofile_pix_by_pix(struct mkonthread *mkp)
if (fabs(array[p]-approx)/array[p] < tolerance)
use_rand_points=0;
+ /* For a check:
+ printf("coord: %g, %g\n", mkp->coord[0], mkp->coord[1]);
+ printf("r_rand: %g (rand: %g, center: %g)\n\n", mkp->r, array[p],
+ approx);
+ */
+
/* Save the peak flux if this is the first pixel: */
if(ispeak) { mkp->peakflux=array[p]; ispeak=0; }
@@ -440,6 +457,10 @@ oneprofile_pix_by_pix(struct mkonthread *mkp)
/* Find the value for this pixel: */
array[p]=profile(mkp);
+ /* For a check:
+ printf("r_center: %g\n", mkp->r);
+ */
+
/* Save the peak flux if this is the first pixel: */
if(ispeak) { mkp->peakflux=array[p]; ispeak=0; }
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 6b8759c..c850428 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -24460,22 +24460,23 @@ vertical).
@deftypefun void gal_box_bound_ellipse_extent (double @code{a}, double
@code{b}, double @code{theta_deg}, double @code{*extent})
Return the maximum extent along each dimension of the given ellipse from
the center of the ellipse. Therefore this is half the extent of the box in
-each dimension. @code{a} is the ellipse major axis, @code{b} is the minor
-axis, @code{theta_deg} is the position angle in degrees. The extent in each
-dimension is in floating point format and stored in @code{extent} which
-must already be allocated before this function.
+each dimension. @code{a} is the ellipse semi-major axis, @code{b} is the
+semi-minor axis, @code{theta_deg} is the position angle in degrees. The
+extent in each dimension is in floating point format and stored in
address@hidden which must already be allocated before this function.
@end deftypefun
@deftypefun void gal_box_bound_ellipse (double @code{a}, double @code{b},
double @code{theta_deg}, long @code{*width})
-Any ellipse can be enclosed into a rectangular box. The purpose of this
-function is to give the height and width of that box assuming the center of
-the ellipse is in the box center. @code{a} is the ellipse major axis,
address@hidden is the minor axis, @code{theta_deg} is the position angle in
-degrees. The @code{width} array will contain the output size in long
-integer type. @code{width[0]}, and @code{width[1]} are the number of pixels
-along the first and second FITS axis. Since the ellipse center is assumed
-to be in the center of the box, all the values in @code{width} will be an
-odd integer.
+Any ellipse can be enclosed into a rectangular box. This function will
+write the height and width of that box where @code{width} points to. It
+assumes the center of the ellipse is located within the central pixel of
+the box. @code{a} is the ellipse semi-major axis length, @code{b} is the
+semi-minor axis, @code{theta_deg} is the position angle in degrees. The
address@hidden array will contain the output size in long integer
+type. @code{width[0]}, and @code{width[1]} are the number of pixels along
+the first and second FITS axis. Since the ellipse center is assumed to be
+in the center of the box, all the values in @code{width} will be an odd
+integer.
@end deftypefun
@deftypefun void gal_box_border_from_center (double @code{center}, size_t
@code{ndim}, long @code{*width}, long @code{*fpixel}, long @code{*lpixel})
diff --git a/lib/box.c b/lib/box.c
index 9d80d37..55f77c6 100644
--- a/lib/box.c
+++ b/lib/box.c
@@ -92,8 +92,8 @@ gal_box_bound_ellipse(double a, double b, double theta_deg,
long *width)
want the final height and width of the box enclosing the
ellipse. So we have to multiply them by two, then take one from
them (for the center). */
- width[0] = 2 * ( (long)extent[0] + 1 ) + 1;
- width[1] = 2 * ( (long)extent[1] + 1 ) + 1;
+ width[0] = 2 * ( (long)extent[0] ) + 1;
+ width[1] = 2 * ( (long)extent[1] ) + 1;
}
diff --git a/lib/gnuastro-internal/kernel-2d.h
b/lib/gnuastro-internal/kernel-2d.h
index ccc2d49..3462467 100644
--- a/lib/gnuastro-internal/kernel-2d.h
+++ b/lib/gnuastro-internal/kernel-2d.h
@@ -34,22 +34,21 @@ along with Gnuastro. If not, see
<http://www.gnu.org/licenses/>.
Make the kernel
---------------
- We'll first make the kernel using MakeProfiles, then crop the outer
- layers that are zero. Put these lines into a script and run it.
+ We'll first make the kernel using MakeProfiles with the following
+ commands. IMPORTANT NOTE: because the kernel is so sharp, random
+ sampling is going to be used for all the pixels (the center won't be
+ used). So it is important to have a large number of random points to
+ make the very slight differences between symmetric parts of the profile
+ even less significant.
- set -o errexit # Stop if a program returns false.
- export GSL_RNG_TYPE=ranlxs2
export GSL_RNG_SEED=1
- astmkprof --kernel=gaussian,2,5 --oversample=1 --envseed -ok2.fits
- astcrop k2.fits --section=2:*-1,2:*-1 --zeroisnotblank \
- --mode=img --output=fwhm2.fits
- rm k2.fits
+ export GSL_RNG_TYPE=ranlxs2
+ astmkprof --kernel=gaussian,2,5 --oversample=1 --envseed
--numrandom=100000
Convert it to C code
--------------------
- We'll use this tiny C program to convert the FITS file into the two
- important C variables:
+ Put the following C program into a file called `kernel.c'.
#include <stdio.h>
#include <stdlib.h>
@@ -60,7 +59,7 @@ along with Gnuastro. If not, see
<http://www.gnu.org/licenses/>.
{
size_t i;
float *arr;
- gal_data_t *img=gal_fits_img_read_to_type("fwhm2.fits", "1",
+ gal_data_t *img=gal_fits_img_read_to_type("kernel.fits", "1",
GAL_TYPE_FLOAT32, -1);
arr=img->array;
@@ -97,33 +96,34 @@ along with Gnuastro. If not, see
<http://www.gnu.org/licenses/>.
-----------------
We can now compile and run that C program and put the outputs in
- `ready.c'. Once its created, copy the contents of `ready.c' after these
- comments.
+ `kernel.c'. Once its created, copy the contents of `kernel-2d.h' after
+ these comments.
- $ astbuildprog -q prepare.c > ready.c
+ $ astbuildprog -q kernel.c > kernel-2d.h
*/
size_t kernel_2d_dsize[2]={11, 11};
-float kernel_2d[121]={0, 0, 0, 0, 0, 2.570758e-08, 0, 0, 0, 0, 0,
+float kernel_2d[121]={0, 0, 0, 0, 0, 2.599797e-08, 0, 0, 0, 0, 0,
+
+0, 0, 3.008479e-08, 6.938075e-07, 4.493532e-06, 8.276223e-06, 4.515019e-06,
6.947793e-07, 3.04628e-08, 0, 0,
-0, 0, 2.981546e-08, 7.249833e-07, 4.468747e-06, 8.409227e-06, 4.554846e-06,
7.034199e-07, 3.002102e-08, 0, 0,
+0, 3.009687e-08, 2.556034e-06, 5.936867e-05, 0.0003808578, 0.0007126221,
0.0003827095, 5.902729e-05, 2.553342e-06, 2.978137e-08, 0,
-0, 3.054498e-08, 2.614154e-06, 5.891601e-05, 0.0003810036, 0.000708165,
0.0003842406, 5.963722e-05, 2.618934e-06, 2.990584e-08, 0,
+0, 7.021852e-07, 5.912285e-05, 0.00137637, 0.008863639, 0.01648383,
0.008855942, 0.001365171, 5.925718e-05, 7.021184e-07, 0,
-0, 7.199899e-07, 5.801019e-05, 0.001365485, 0.009023659, 0.01638159,
0.008892864, 0.001345278, 5.920425e-05, 6.984741e-07, 0,
+0, 4.490787e-06, 0.0003826718, 0.008857355, 0.05742518, 0.1062628, 0.05727194,
0.008880079, 0.0003826067, 4.478989e-06, 0,
-0, 4.584869e-06, 0.0003830431, 0.008917402, 0.05743425, 0.1061489, 0.05746412,
0.008902563, 0.0003849257, 4.448404e-06, 0,
+2.595735e-08, 8.31301e-06, 0.0007113572, 0.01640853, 0.1061298, 0.1971036,
0.1062611, 0.01647962, 0.000708363, 8.379878e-06, 2.593496e-08,
-2.572769e-08, 8.414041e-06, 0.0007008284, 0.0164456, 0.1055995, 0.19753,
0.1061855, 0.01653461, 0.0007141303, 8.41643e-06, 2.550312e-08,
+0, 4.516684e-06, 0.0003846966, 0.008860709, 0.05739478, 0.1062216, 0.05725683,
0.00881713, 0.000383981, 4.473017e-06, 0,
-0, 4.582525e-06, 0.0003775396, 0.00898499, 0.05741534, 0.1062144, 0.05700329,
0.008838926, 0.0003822096, 4.543726e-06, 0,
+0, 6.950547e-07, 5.920586e-05, 0.00137483, 0.00887785, 0.0164709, 0.008855232,
0.001372743, 5.939038e-05, 7.016624e-07, 0,
-0, 6.883925e-07, 6.09077e-05, 0.001339333, 0.008817007, 0.01636454,
0.008995386, 0.001407854, 6.004799e-05, 7.203602e-07, 0,
+0, 3.006322e-08, 2.587011e-06, 5.92911e-05, 0.0003843824, 0.0007118155,
0.000386519, 5.974654e-05, 2.585581e-06, 3.048036e-08, 0,
-0, 3.095966e-08, 2.575403e-06, 5.89859e-05, 0.0003804447, 0.0007091904,
0.0003810006, 5.903253e-05, 2.575202e-06, 2.934356e-08, 0,
+0, 0, 3.041056e-08, 7.05225e-07, 4.497418e-06, 8.388542e-06, 4.478833e-06,
7.018358e-07, 2.995504e-08, 0, 0,
-0, 0, 3.040937e-08, 7.018197e-07, 4.543086e-06, 8.296753e-06, 4.434901e-06,
6.659026e-07, 3.066215e-08, 0, 0,
+0, 0, 0, 0, 0, 2.567377e-08, 0, 0, 0, 0, 0};
-0, 0, 0, 0, 0, 2.603901e-08, 0, 0, 0, 0, 0};
#endif
[Prev in Thread] |
Current Thread |
[Next in Thread] |
- [gnuastro-commits] master b4c2b5a: gal_box_bound_ellipse return width isn't extra,
Mohammad Akhlaghi <=