[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[gnuastro-commits] master 0dec9ac 110/125: Allow for floating point erro
From: |
Mohammad Akhlaghi |
Subject: |
[gnuastro-commits] master 0dec9ac 110/125: Allow for floating point errors in Crop comparisons |
Date: |
Sun, 23 Apr 2017 22:36:49 -0400 (EDT) |
branch: master
commit 0dec9ac3b0af3ed318570c33d927cb9941044115
Author: Mohammad Akhlaghi <address@hidden>
Commit: Mohammad Akhlaghi <address@hidden>
Allow for floating point errors in Crop comparisons
When the input image is not aligned with the WCS or the pixel scales along
each dimension are different, in WCS mode, Crop will complain and abort
(because it currently uses these assumptions). Until now this check was
done exactly, so even a small difference of, for example 1e-10 (which is
clearly due to floating point errors) would cause an error. So with this
commit, very small relative differences (less than 1e-6) are acceptable and
Crop won't complain.
In this commit, short explanation that was also added to the `--polygon'
option of Crop so users can use SAO ds9 (a FITS viewer) for defining the
polygon instead of having to type the values manually.
Both these suggestions were kindly made by Guillaume Mahler, so his name
has been added to the `THANKS' file.
---
THANKS | 1 +
bin/crop/crop.c | 4 ++++
bin/crop/wcsmode.c | 33 ++++++++++++++++++---------------
doc/gnuastro.texi | 43 +++++++++++++++++++++++++++++++++++++------
4 files changed, 60 insertions(+), 21 deletions(-)
diff --git a/THANKS b/THANKS
index b1d13c7..a28ae7c 100644
--- a/THANKS
+++ b/THANKS
@@ -23,6 +23,7 @@ support in Gnuastro. The list is ordered alphabetically.
Brandon Invergo address@hidden
Lee Kelvin address@hidden
Mohammad-Reza Khellat address@hidden
+ Guillaume Mahler address@hidden
Alan Lefor address@hidden
Francesco Montanari address@hidden
William Pence address@hidden
diff --git a/bin/crop/crop.c b/bin/crop/crop.c
index 12a708e..2de990d 100644
--- a/bin/crop/crop.c
+++ b/bin/crop/crop.c
@@ -308,6 +308,10 @@ wcsmodecrop(void *inparam)
}
while ( ++(crp->in_ind) < p->numin );
+ /* Correct in_ind. The loop above went until `in_ind' is one more
+ than the index for the last input image (that is how it exited the
+ loop). But `crp->in_ind' is needed later, so correct it here. */
+ --crp->in_ind;
/* Check the final output: */
if(crp->numimg)
diff --git a/bin/crop/wcsmode.c b/bin/crop/wcsmode.c
index 3173ccf..77f3684 100644
--- a/bin/crop/wcsmode.c
+++ b/bin/crop/wcsmode.c
@@ -56,15 +56,18 @@ wcs_check_prepare(struct cropparams *p, struct inputimgs
*img)
double imgcrd[8], phi[4], theta[4], pixcrd[8];
- /* Check if the image is aligned correctly and fits the
- resolution of other images. */
- if(wcs->pc[1]!=0.0f || wcs->pc[2]!=0.0f)
+ /* Check if the image is aligned with the WCS coordinates. Note that
+ because of small floating point errors, some programs might still keep
+ very small values in the off-diagonal matrix elements. */
+ if( fabs(wcs->pc[1]/wcs->pc[3])>1e-6 || fabs(wcs->pc[2]/wcs->pc[3])>1e-6 )
error(EXIT_FAILURE, 0, "%s: HDU %s: is not aligned to the "
"celestial coordinates. The first FITS axis should be "
"along the Right Ascension and the second FITS axis "
- "should be along the declination. You should rotate "
- "(interpolate) the images with other software",
- img->name, p->cp.hdu);
+ "should be along the declination.\n\n"
+ "Gnuastro's Warp program can align it with the following "
+ "command:\n\n"
+ " $ astwarp %s --hdu=%s --align\n",
+ img->name, p->cp.hdu, img->name, p->cp.hdu);
if(wcs->pc[0]>0)
error(EXIT_FAILURE, 0, "%s: HDU %s: An increase in the first "
"FITS axis pixel coordinates should be a decrese in the "
@@ -78,19 +81,20 @@ wcs_check_prepare(struct cropparams *p, struct inputimgs
*img)
img->name, p->cp.hdu);
- /* Since we are dealing with very accurate values, a multiplication
- by -1 might cause a floating point error. So we have to account
- for the floating point error. */
- if(-1.0f*wcs->pc[0]<wcs->pc[3]-1e-15 || -1.0f*wcs->pc[0]>wcs->pc[3]+1e-15)
+ /* Check the image resolution (if the pixels are actually a square), then
+ compare the resolution with the other input images. Due to floating
+ point errors, some very small differences might exist in the pixel
+ scale, so break out with an error only if the pixel scales are more
+ different than 1e-6. */
+ pixscale=gal_wcs_pixel_scale_deg(wcs);
+ if( fabs(pixscale[0]-pixscale[1])/pixscale[0] > 1e-6 )
error(EXIT_FAILURE, 0, "%s: HDU %s: The pixel scale along "
"the two image axises is not the same. The first axis "
- "is %f arcseconds/pixel, while the second is %f",
- img->name, p->cp.hdu, 3600*-1.0f*wcs->pc[0],
- 3600*wcs->pc[3]);
+ "is %.15g deg/pixel, while the second is %.15g",
+ img->name, p->cp.hdu, pixscale[1], pixscale[0]);
if(p->res==0.0f)
{
/* Set the resolution of the image. */
- pixscale=gal_wcs_pixel_scale_deg(wcs);
p->res=pixscale[0];
free(pixscale);
@@ -114,7 +118,6 @@ wcs_check_prepare(struct cropparams *p, struct inputimgs
*img)
}
else
{
- pixscale=gal_wcs_pixel_scale_deg(wcs);
if(p->res!=pixscale[0])
error(EXIT_FAILURE, 0, "%s: HDU %s: The resolution of "
"this image is %f arcseconds/pixel while the "
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 0a8c8de..9717e7f 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -7853,14 +7853,14 @@ Crop box parameters:
@item -x FLT
@itemx --xc=FLT
-X axis (first image axis, horizontal when viewed in SAO DS9) position of
+X axis (first image axis, horizontal when viewed in SAO ds9) position of
crop box center. Along with @option{--yc}, this only produces one crop in
each run. The width of the square crop (in pixels) is the value to
@option{--iwidth}.
@item -y FLT
@itemx --yc=FLT
-Y axis (second image axis, vertical when viewed in SAO DS9) position of
+Y axis (second image axis, vertical when viewed in SAO ds9) position of
crop box center. See @option{--xc}.
@item -r FLT
@@ -7929,10 +7929,41 @@ have set the operating mode to WCS mode in your nearest
configuration file
images/tiles and Crop will stitch them to produce this cropped region:
@example
-$ astcrop --mode=wcs desired-filter-image(s).fits \
- --polygon="53.187414,-27.779152 : 53.159507,-27.759633 : \
+$ astcrop --mode=wcs desired-filter-image(s).fits \
+ --polygon="53.187414,-27.779152 : 53.159507,-27.759633 : \
53.134517,-27.787144 : 53.161906,-27.807208"
@end example
address@hidden SAO ds9
+In other cases, you have an image and want to define the polygon yourself
+(it isn't already published like the example above). As the number of
+vertices increases, checking the vertice coordinates on a FITS viewer (for
+example SAO ds9) and typing them in one by one can be very tedious and
+prone to typo errors.
+
+You can take the following steps to avoid the frustration and possible
+typos: Open the image with ds9 and activate its ``region'' mode with
address@hidden@click{}Region}. Then define the region as a polygon
+with @address@hidden@click{}Polygon}. Click on the
+approximate center of the region you want and a small square will
+appear. By clicking on the vertices of the square you can shrink or expand
+it, clicking and dragging anywhere on the edges will enable you to define a
+new vertice. After the region has been nicely defined, save it as a file
+with @address@hidden Regions}. You can then select the
+name and address of the output file, keep the format as @command{REG} and
+press ``OK''. In the next window, keep format as ``ds9'' and ``Coordinate
+System'' as ``fk5''. A plain text file (let's call it @file{ds9.reg}) is
+now created.
+
+You can now convert this plain text file to Crop's polygon format with this
+command (when typing on the command-line, ignore the address@hidden'' at the
end
+of the first and second lines along with the extra spaces, these are only
+for nice printing):
+
address@hidden
+$ v=$(awk 'NR==4' ds9.reg | sed -e's/polygon(//' \
+ -e's/\([^,]*,[^,]*\),/\1:/g' -e's/)//' )
+$ astcrop --mode=wcs image.fits --polygon=$v
address@hidden example
@item --outpolygon
Keep all the regions outside the polygon and mask the inner ones with blank
@@ -10601,7 +10632,7 @@ the examples above as a demonstration.
@cindex FITS standard
Based on the FITS standard, integer values are assigned to the center of a
pixel and the coordinate [1.0, 1.0] is the center of the first pixel
-(bottom left of the image when viewed in SAO DS9). So the coordinate center
+(bottom left of the image when viewed in SAO ds9). So the coordinate center
[0.0, 0.0] is half a pixel away (in each axis) from the bottom left vertice
of the first pixel. The resampling that is done in Warp (see
@ref{Resampling}) is done on the coordinate axises and thus directly
@@ -10716,7 +10747,7 @@ ignored.
@item -c
@itemx --centeroncorer
Put the center of coordinates on the corner of the first (bottom-left when
-viewed in SAO DS9) pixel. This option is applied after the final warping
+viewed in SAO ds9) pixel. This option is applied after the final warping
matrix has been finalized: either through modular warpings or the raw
matrix. See the explanation above for coordinates in the FITS standard to
better understand this option and when it should be used.
- [gnuastro-commits] master 7e888ce 073/125: First implementation of MakeProfiles with gal_data_t, (continued)
- [gnuastro-commits] master 7e888ce 073/125: First implementation of MakeProfiles with gal_data_t, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master e8d9751 087/125: MakeNoise now uses the new gal_data_t infrastructure, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master e978d44 068/125: Numeric option sanity checks, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master e74eb25 093/125: Changed name of Header program to Fits, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master e353932 109/125: NaN magnitude will only blank the profile's area, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 96ab804 098/125: Working on tiles implemented in tile checking, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 31b3b76 099/125: Spatial convolution now using new tessellation, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master c0dd6db 102/125: The noedgecorrection option in spatial domain convolution, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 0dfdc34 069/125: Option values directly set in library, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 95ba655 100/125: Spatial convolution with new tessellation, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 0dec9ac 110/125: Allow for floating point errors in Crop comparisons,
Mohammad Akhlaghi <=
- [gnuastro-commits] master db4c16d 076/125: Pulled in all the recent work in master, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 64f8db5 080/125: ConvertType working except for text inputs and outputs, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 4d853b4 112/125: Old neighbors.h header removed, used new in MakeProfiles, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 9959e45 004/125: Arithmetic using new data structure, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 0a2fd0e 089/125: New implementations for histogram and CFP, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 213d1c3 106/125: gal_data_copy_* functions can now work on tiles, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 4f6a8c4 095/125: Fits prints HDU info of FITS file with no options, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 63043e5 105/125: Statistical functions can accept tiled input, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master 4839a11 104/125: Preparations for tile interpolation, Mohammad Akhlaghi, 2017/04/23
- [gnuastro-commits] master d6aa57f 081/125: Text library can write images also, besides tables, Mohammad Akhlaghi, 2017/04/23