[gnuastro-commits] master 24ccb25 2/2: Crop: Minor corrections in concave-sort and its usage in Crop

[Mohammad Akhlaghi]

branch: master
commit 24ccb251cdb50c72dcf7c470cdf4a56040cfa2bd
Author: Mohammad Akhlaghi <address@hidden>
Commit: Mohammad Akhlaghi <address@hidden>

    Crop: Minor corrections in concave-sort and its usage in Crop
    
    In the previous commit, Sachin implemented the concave sort function as a
    library function in Gnuastro and also called it within Crop. Only one issue
    was found regarding the check for the warning of sorted concave
    polygon. The check had to be done after we sort the polygon vertices in a
    final array.
    
    With this commit, that was fixed and some stylistic issues were also
    implemented to make the code blend in better with the rest of Gnuastro's
    source. The documentation of `gal_polygon_vertices_sort' was also improved
    in the book and an example was added.
    
    This completes task #13658.
---
 NEWS               |  6 +++++
 bin/crop/onecrop.c | 36 +++++++++++++++++++-----------
 doc/gnuastro.texi  | 28 ++++++++++++++++++-----
 lib/polygon.c      | 65 +++++++++++++++++++++++-------------------------------
 4 files changed, 79 insertions(+), 56 deletions(-)

diff --git a/NEWS b/NEWS
index 70d50c7..ca9031b 100644
--- a/NEWS
+++ b/NEWS
@@ -62,6 +62,12 @@ See the end of the file for license conditions.
    - gal_polygon_is_inside: if point is inside polygon (convex or concave).
    - gal_polygon_is_counterclockwise: check if polygon is counter-clockwise.
    - gal_polygon_to_counterclockwise: convert to counter-clockwise if it isn't.
+   - gal_polygon_vertices_sort: un-ordered vertices to concave/convext 
polygons.
+   - gal_units_extract_decimal: Extract numbers from strings like "A:B:C".
+   - gal_units_ra_to_degree: Convert RA (HH:MM:SS) to degrees.
+   - gal_units_dec_to_degree: Convert Dec (DD:MM:SS) to degrees.
+   - gal_units_degree_to_ra: Convert degrees to RA (DD:MM:SS).
+   - gal_units_degree_to_dec: Convert degrees to Dec (DD:MM:SS).
 
 ** Removed features
 
diff --git a/bin/crop/onecrop.c b/bin/crop/onecrop.c
index 602b5cf..24343a2 100644
--- a/bin/crop/onecrop.c
+++ b/bin/crop/onecrop.c
@@ -363,37 +363,47 @@ polygonmask(struct onecropparams *crp, void *array, long 
*fpixel_i,
     error(EXIT_FAILURE, errno, "%s: allocating %zu bytes for ordinds",
           __func__, nvertices*sizeof *ordinds);
 
-  /* If the user wants to sort the edges, do it, if not, make sure its in
+  /* If the user wants to sort the edges, do it. If not, make sure its in
      counter-clockwise orientation. */
   if(crp->p->polygonsort)
-    {
-      gal_polygon_vertices_sort(crp->ipolygon, crp->p->nvertices, ordinds);
-
-      if( !gal_polygon_is_convex(crp->ipolygon, crp->p->nvertices))
-        error(0, 0, "%s: Warning: The sorted concave polygon might "
-              "not be in the same desired format as required",
-              __func__);
-    }
+    gal_polygon_vertices_sort(crp->ipolygon, nvertices, ordinds);
   else
     {
-      for(i=0;i<crp->p->nvertices;++i) ordinds[i]=i;
-      gal_polygon_to_counterclockwise(crp->ipolygon, crp->p->nvertices);
+      /* Keep the original order of the vertices, just make it
+         counter-clockwise. */
+      for(i=0;i<nvertices;++i) ordinds[i]=i;
+      gal_polygon_to_counterclockwise(crp->ipolygon, nvertices);
     }
 
   /* Fill the final polygon vertice positions within `ipolygon' and also
      the fpixel_i coordinates from all the vertices to bring them into the
      crop image coordinates. */
-  for(i=0;i<crp->p->nvertices;++i)
+  for(i=0;i<nvertices;++i)
     {
       ipolygon[i*2  ] = crp->ipolygon[ordinds[i]*2]   - fpixel_i[0];
       ipolygon[i*2+1] = crp->ipolygon[ordinds[i]*2+1] - fpixel_i[1];
     }
 
+  /* Print a warning if we done the sorting, _and_ the sorted polygon is
+     concave, _and_ the user hasn't activated the quiet mode. Note that we
+     could'n do this immediately after calling `gal_polygon_vertices_sort'
+     because that function doesn't touch the actual vertice values, it only
+     fills `ordinds'. */
+  if(crp->p->polygonsort
+     && !crp->p->cp.quiet
+     && !gal_polygon_is_convex(ipolygon, nvertices) )
+    error(0, 0, "%s: warning: the given vertices belong to a concave "
+          "polygon, but there is no unique solution to the sorting of "
+          "concave polygons. Please check the cropped image, if it doesn't "
+          "correspond to your desired polygon, sort the vertices yourself "
+          "in counter-clockwise order _and_ don't use the `--polygonsort' "
+          "option", __func__);
+
   /* Set the function for checking if a point is inside the polygon. For
      concave polygons the process is more complex and thus
      slower. Therefore when the polygon is convex, its better to use the
      simpler/faster function. */
-  isinside = ( gal_polygon_is_convex(ipolygon, crp->p->nvertices)
+  isinside = ( gal_polygon_is_convex(ipolygon, nvertices)
                ? gal_polygon_is_inside_convex
                : gal_polygon_is_inside );
 
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index db7b569..6e7d974 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -23303,12 +23303,28 @@ This function uses the 
@url{https://en.wikipedia.org/wiki/Sutherland%E2%80%93Hod
 Note that the vertices of both polygons have to be sorted in an 
anti-clock-wise manner.
 @end deftypefun
 
-@deftypefun void gal_polygon_vertices_sort(double *in, size_t n, size_t 
*ordinds)
-The sorting function which takes points as inputs and makes a sorted array 
(@code(ordinds}) of indexes.
-In this algorithm, we find the rightmost ans leftmost points(based on their 
x-coordinate) and use the diagonal vector between those points to group the 
points in arrays based on their position with respect to this vector.
-Now, for anticlockwise sorting, all the points below the vector are sorted in 
ascending of their x-coordinates and points above the vector are sorted in 
decreasing order using @code{qsort}.
-Finally, both these arrays are merged together to get the final sorted array 
of points, from which the points are indexed into the@code{ordinds} using 
linear search.
-Remember as there is no unique way to sort concave polygons the final sorted 
output may not be the one desired.
+@deftypefun void gal_polygon_vertices_sort (double @code{*vertices}, size_t 
@code{n}, size_t @code{*ordinds})
+Sort the indexs of the un-ordered @code{vertices} array to a counter-clockwise 
polygon in the already allocated space of @code{ordinds}.
+It is assumed that there are @code{n} vertices, and thus that @code{vertices} 
contains @code{2*n} elements where the two coordinates of the first vertice 
occupy the first two elements of the array and so on.
+
+The polygon can be both concave and convex (see the start of this section).
+However, note that for concave polygons there is no unique sort from an 
un-ordered set of vertices.
+So after this function you may want to use @code{gal_polygon_is_convex} and 
print a warning to check the output if the polygon was concave.
+
+Note that the contents of the @code{vertices} array are left untouched by this 
function.
+If you want to write the ordered vertice coordinates in another array with the 
same size, you can use a loop like this:
+
+@example
+for(i=0;i<n;++i)
+@{
+  ordered[i*2  ] = vertices[ ordinds[i]*2    ];
+  ordered[i*2+1] = vertices[ ordinds[i]*2 + 1];
+@}
+@end example
+
+In this algorithm, we find the rightmost and leftmost points (based on their 
x-coordinate) and use the diagonal vector between those points to group the 
points in arrays based on their position with respect to this vector.
+For anticlockwise sorting, all the points below the vector are sorted by their 
ascending x-coordinates and points above the vector are sorted in decreasing 
order using @code{qsort}.
+Finally, both these arrays are merged together to get the final sorted array 
of points, from which the points are indexed into the @code{ordinds} using 
linear search.
 @end deftypefun
 
 
diff --git a/lib/polygon.c b/lib/polygon.c
index 7729a60..639df99 100644
--- a/lib/polygon.c
+++ b/lib/polygon.c
@@ -807,29 +807,30 @@ polygon_make_arr(double *in, size_t n, size_t *A_size, 
size_t *B_size,
 
 
 
-/* The comparator functions for qsort. CompareA arranges
-    the array in ascending order according to their
-    x-coordinate and CompareB arranges in descending
-    order according to their x-coordintes.
-    */
+/* The comparator functions for qsort. CompareA arranges the array in
+   ascending order according to their x-coordinate. */
 static int
 polygon_compareA(const void *a, const void *b)
 {
   struct point *p1 = (struct point *)a, *p2 = (struct point *)b;
-
-  return (p1->x==p2->x) ? 0:( (p1->x<p2->x) ? -1:1 );
+  return ( p1->x==p2->x
+           ? 0
+           : (p1->x<p2->x ? -1 : 1) );
 }
 
 
 
 
 
+/* The comparator functions for qsort. CompareB arranges in descending
+   order according to their x-coordintes. */
 static int
 polygon_compareB(const void *a, const void *b)
 {
   struct point *p1 = (struct point *)a, *p2 = (struct point *)b;
-
-  return (p1->x==p2->x) ? 0:( (p1->x<p2->x) ? 1:-1 );
+  return ( p1->x==p2->x
+           ? 0
+           : (p1->x<p2->x ? 1 : -1) );
 }
 
 
@@ -840,7 +841,7 @@ polygon_compareB(const void *a, const void *b)
     together to the output array. Hence it is the main function
     which should be called when using concave sort. */
 void
-gal_polygon_vertices_sort(double *in, size_t n, size_t *ordinds)
+gal_polygon_vertices_sort(double *vertices, size_t n, size_t *ordinds)
 {
   size_t i, j, A_size = 0, B_size = 0;
   struct point A[GAL_POLYGON_MAX_CORNERS];
@@ -848,7 +849,7 @@ gal_polygon_vertices_sort(double *in, size_t n, size_t 
*ordinds)
   struct point sorted[GAL_POLYGON_MAX_CORNERS];
   struct point tordinds[GAL_POLYGON_MAX_CORNERS];
 
-
+  /* Sanity check. */
   if(n>GAL_POLYGON_MAX_CORNERS)
     error(EXIT_FAILURE, 0, "%s: most probably a bug! The number of corners "
           "is more than %d. This is an internal value and cannot be set from "
@@ -856,42 +857,32 @@ gal_polygon_vertices_sort(double *in, size_t n, size_t 
*ordinds)
           "value. Please contact us at %s so we can solve this problem",
           __func__, GAL_POLYGON_MAX_CORNERS, PACKAGE_BUGREPORT);
 
+  /* Make arrays A and B and store the vertices in them. Currently points
+     are stored in based on their position from the diagonal vector. */
+  polygon_make_arr(vertices, n, &A_size, &B_size, A, B);
 
-  /* Make arrays A and B and store the vertices in them.
-     Currently points are stored in based on their position
-     from the diagonal vector.
-    */
-
-
-  polygon_make_arr(in, n, &A_size, &B_size, A, B);
-
-  /* Now, we put the contents of A and B in the temporary array.
-     Firstly, we put the contents of A and then save the last index
-     of A(stored in i) and continue from that index while copying
-     from B(using j).
-  */
+  /* Now, we put the contents of A and B in the temporary array. Firstly,
+     we put the contents of A and then save the last index of A (stored in
+     i) and continue from that index while copying from B(using j). */
   for(i=0; i<A_size; i++) tordinds[i]=A[i];
   for(j=0; j<B_size; j++) tordinds[i++]=B[j];
 
-  /* Now sort the arrays A and B w.r.t their x axis,
-     sorting A in ascending order and B in descending order. */
+  /* Now sort the arrays A and B w.r.t their x axis, sorting A in ascending
+     order and B in descending order. */
   qsort(A, A_size, sizeof(struct point), polygon_compareA);
   qsort(B, B_size, sizeof(struct point), polygon_compareB);
 
-/*Finally, we put the contents of A and B in a final sorted array.*/
+  /*Finally, we put the contents of A and B in a final sorted array.*/
   for(i=0; i<A_size; i++) sorted[i]=A[i];
   for(j=0; j<B_size; j++) sorted[i++]=B[j];
 
-
-/* The temporary array is now used to find the location of points
-   stored in sorted array and assign index in ordinds accordingly.*/
+  /* The temporary array is now used to find the location of points stored
+     in sorted array and assign index in ordinds accordingly.*/
   for(i=0; i<n; i++)
     for(j=0; j<n; j++)
-        if( (tordinds[i].x == sorted[j].x) &&
-            (tordinds[i].y == sorted[j].y) )
-          {
-            ordinds[j]=i;
-            break;
-          }
-
+      if( tordinds[i].x == sorted[j].x && tordinds[i].y == sorted[j].y )
+        {
+          ordinds[j]=i;
+          break;
+        }
 }