taskGNU Astronomy Utilities - Tasks: task #13658, Work on concave polygons too

 
 

You are not allowed to post comments on this tracker with your current authentication level.

task #13658: Work on concave polygons too

Submitter:  Mohammad Akhlaghi <makhlaghi>
Submitted:  Thu 25 Jun 2015 02:44:12 AM UTC
   
 
Should Start On:  Wed 24 Jun 2015 03:00:00 PM UTC Should be Finished on:  Wed 24 Jun 2015 03:00:00 PM UTC
Category:  Crop Priority:  5 - Normal
Item Group:  New feature Status:  Done
Privacy:  Public Assigned to:  sks_15
Percent Complete:  100% Open/Closed:  Closed
Effort:  0.00

Jump to the original submission

Thu 09 Apr 2020 04:57:33 PM UTC, comment #48: 

Thanks a lot Sachin, it has been merged and pushed to the main branch :-D.

Just some tips for the future (which I have implemented in the next commit):

  • In the documentation, when you define a function, put name of the argument variables in a `@code{}'.


  • When printing a warning in the programs, check if the `--quiet' option has been called by the user. This is part of the common options, and in the code, its in the structure `gal_options_common_params' (which is defined in lib/gnuastro-internal/options.h. In all the programs this structure is in `p->cp', so you can check the `--quiet' option with `p->cp.quiet'.


  • The precedence of operators like `&&' and `?:' are very low in the C language, so there is no need to use parenthesis when they are placed between comparison operators like `<' or `=='. Parenthesis (like brackets) should only be used when the precedence is different from the C standard. An extra parenthesis will just make the code-reader wonder why you used them? For example, they are going to ask themselvs: Did you intend to do something else and forgot? is it a bug? The Parenthesis I have added at the start and and end of the return statements of `polygon_compareA' and `polygon_compareB' are to help in automatic formatting of the code when it is broken up in multiple lines. Usually as soon as the `?:' statement gets a little complex, its much more easier to read if you break it into multiple lines like these two functions.


  • Try to set descriptive names for your arguments and avoid using something like `in' or `input', because they are very vague. For example in the case of `gal_polygon_vertices_sort', I have changed the first argument to be `vertices' instead of `in'. Its a form of documentation, and also makes it easier to search for the instances of the variable in the code. It also helps in actually writing the description of the function in the book ;-).


Mohammad Akhlaghi <makhlaghi>
Group administrator
Thu 09 Apr 2020 07:14:11 AM UTC, comment #47: 

Done here:-)

Sachin Kumar Singh <sks_15>
Thu 09 Apr 2020 03:10:12 AM UTC, comment #46: 

Thanks Sachin. I just made another merge for task #15567, it would be cleaner if you re-base your commit over the updated master branch, then I'll merge it ;-). You can also make the final touches ;-). I won't make any merge or commits to the master branch till I get your commit (if you do it by the end of Thursday) ;-).

Mohammad Akhlaghi <makhlaghi>
Group administrator
Wed 08 Apr 2020 07:16:47 PM UTC, comment #45: 

Done here. Also, sorry I forgot about the usage of conditional operators for if else in compare functions before rebasing and pushing. But probably you will be doing some more changes in a new commit and hence I didn't rebase again. Can you please fix them when you do the final changes? Thank you;-)

Sachin Kumar Singh <sks_15>
Wed 08 Apr 2020 02:01:40 PM UTC, comment #44: 

Great! Thanks Sachin. A few small points remain, I look forward to merging it afterwards:

  • The commit isn't on the most recent commit of the master branch. You can fetch the most recent commit(s) on the master branch and then rebase over them.


  • The static structures and functions that are only defined in `polygon.c' (and not `polygon.h') are only relevant during the building of the library. They aren't accessible for a user of Gnuastro once it is installed. So they don't need documentation in the book. They should only be documented in the actual source as comments.


  • In the documentation of `gal_polygon_vertices_sort', it would be really helpful if you explain the methodology of your algorithm a little. In particular note that for concave polygons there is no unique sorting and that people should beware (possibly using the polygon-type function to make sure if it is a convex or concave polygon afterwards). Also, do this when you call this function in the Crop program's `onecrop.c': put a check after the sorting and if the polygon is convex, print a warning, with `error(0,0,"WARNING: ...."')  and let the user know that the sorted convex polygon may not be what they wanted.


  • You can set `struct point' to be static in `polygon.c'.


  • I noticed multiple return statements in `polygon_compareA' and `polygon_compareB' that can be made more readable with the `?:' structure ;-).


  • Finally, two extra curly-brace blocks (after the `for's) can be removed in the final set of loops in `gal_polygon_vertices_sort'. Extra curly braces can be buggy and will always confuse the reader on why you used them? For example was there something you wanted to change but forgot?
Mohammad Akhlaghi <makhlaghi>
Group administrator
Wed 08 Apr 2020 12:22:56 PM UTC, comment #43: 

All done:-) Documentation is done and rebased here. One minor point, when I tried to declare `struct point p` in polyon_leftmost_point` and used p={..}, there was the error that an expression has to be provided(it came in the beginning and that's why I declared it in the middle to start with). Hence I haven't change it. Other than that it's all good:-)

Sachin Kumar Singh <sks_15>
Tue 07 Apr 2020 10:15:46 PM UTC, comment #42: 

One more minor stylistic point:

When you have a simple test, its much elegant/readable/less-buggy to use the `?:' structure, rather than repeating something.

For the return value of `polygon_leftof_vector', you have such a statement. So its better to write it like this:


return test ? (test>0 ? 1 : -1) : 0;


Also, please define variables at the top of functions, not in the middle of its operations (as discussed in Comment 55 of Task 15567.

Mohammad Akhlaghi <makhlaghi>
Group administrator
Tue 07 Apr 2020 09:56:56 PM UTC, comment #41: 

Great! Thanks a lot Sachin, you have the green flag :-). Just some minor points are listed below to hopefully implement with the rebase:

  • When including headers, and order doesn't matter for the algorithm/compilation, please sort them by length. So in `polygon.c', include the `float.h' before `stdlib.h'. Ideally, when the lengths are also equal, it would also be more pleasing to sort them alphabetically ;-). The basic C library headeres are usually independent of each other.


  • It would help in readability if you don't include extra curly braces (`{' or `}'). For example in `polygon_rightmost_point', the `if' statement is right after the `for' statement and there is no other operation. So you can just write it like below. Also note that as discussed above, because the `max_index' line is independent of (doesn't affect) the `tmp_max' line, and its shorter, its better for it to be above the `tmp_max' line ;-). I also see a few such cases in `gal_polygon_vertices_sort' (the two `i' and `j' loops, for the latter, you can just write it like this to be more readable: `tordinds[i++]=B[j];').



  for(i=0; i<n; i++)
    if(tmp_max < in[i*2])
      {
        max_index = i;
        tmp_max = in[i*2];
      }


  • Some static functions don't have a `polygon_' prefix (like `compareA' and `compareB').
Mohammad Akhlaghi <makhlaghi>
Group administrator
Mon 06 Apr 2020 10:37:24 PM UTC, comment #40: 

All things are done according to the requirements. Also the `gal_polygon_sort_vertex` function, until now, returned the array of output points rather than indexes(as ordinds). The latest commit fixes that and puts point structure in polygon.c.(I also toyed with the idea of a map to sort, but it turned out to be rather redundant and overkilling so was dropped altogether).

For storing indexes in ordinds, I used a linear search algorithm giving a worst case of O(n^2). A binary search (O(log2(n))seems to be a better alternative but as the size of array won't be very large(GAL_POLYGON_MAX_CORNERS == 50) I think a linear search won't be that bad:-).

Documentation and a slight error(a silly whitespace) will be fixed and rebased as soon as the final commit gets a green flag:-)

Sachin Kumar Singh <sks_15>
Sun 05 Apr 2020 09:08:54 PM UTC, comment #39: 

Thanks Sachin, it generally looks good, can you please rebase it so I can implement it in the Crop program, and close this task?

Or maybe you can also implement it in Crop yourself? Effectively, just replacing the call to `gal_polygon_vertices_sort_convex' with `gal_polygon_vertices_sort'

I just saw this warning when compiling the branch:


In file included from /path/to/gnuastro/source/bin/crop/onecrop.c:37:
/path/to/gnuastro/source/lib/gnuastro/polygon.h:64:1: warning: useless storage class specifier in empty declaration
   64 | };
      | ^


Can you please fix it? Generally, you should't put `static' things in a header! You can put it directly in `polygon.c', it just has to be above the first instance of its usage ;-).

About Github, you have your own fork there and there is no problem, but Github is not very liked in the free software community because its not free software itself (ironic, isn't it?)! GNU has an Ethical Repository Criteria Evaluations which you can see for more information. For a more modern visual experience, Gitlab is far better than Github.

Generally, like many other GNU packages, Gnuastro's main Git repository is on Savannah's CGIT interface.

Mohammad Akhlaghi <makhlaghi>
Group administrator
Sat 04 Apr 2020 04:56:58 AM UTC, comment #38: 

Yes, I did get this error and proposed for the same solutions in comment #36 which is now fixed here. Also can you please push all the recent changes of savannah master in GitHub master branch too. It helps in maintaining consistency.

Sachin Kumar Singh <sks_15>
Sat 04 Apr 2020 03:20:20 AM UTC, comment #37: 

I tried the build, but it crashed when building `onecrop.c':


gcc -DHAVE_CONFIG_H -I. -I/path/to/gnuastro/source/bin/crop -I../..  -I\/path/to/gnuastro/source/bootstrapped/lib -I\/path/to/gnuastro/source/lib   -Wall -g -O0  -pthread -MT onecrop.o -MD -MP -MF $depbase.Tpo -c -o onecrop.o /path/to/gnuastro/source/bin/crop/onecrop.c &&\
mv -f $depbase.Tpo $depbase.Po
In file included from /path/to/gnuastro/source/bin/crop/onecrop.c:37:
/path/to/gnuastro/source/lib/gnuastro/polygon.h:64:1: warning: useless storage class specifier in empty declaration
   64 | };
      | ^
In file included from /path/to/gnuastro/source/bin/crop/onecrop.c:37:
/path/to/gnuastro/source/lib/gnuastro/polygon.h:102:1: warning: ‘polygon_leftmost_point’ declared ‘static’ but never defined [-Wunused-function]
  102 | polygon_leftmost_point(double *in, size_t n);
      | ^~~~~~~~~~~~~~~~~~~~~~
/path/to/gnuastro/source/lib/gnuastro/polygon.h:105:1: warning: ‘polygon_rightmost_point’ declared ‘static’ but never defined [-Wunused-function]
  105 | polygon_rightmost_point(double *in, size_t n);
      | ^~~~~~~~~~~~~~~~~~~~~~~
/path/to/gnuastro/source/lib/gnuastro/polygon.h:108:1: warning: ‘polygon_leftof_vector’ declared ‘static’ but never defined [-Wunused-function]
  108 | polygon_leftof_vector(double *in, size_t n, double x, double y);
      | ^~~~~~~~~~~~~~~~~~~~~
/path/to/gnuastro/source/lib/gnuastro/polygon.h:111:1: warning: ‘polygon_make_arr’ declared ‘static’ but never defined [-Wunused-function]
  111 | polygon_make_arr(double *in, size_t n, size_t *A_size, size_t *B_size,
      | ^~~~~~~~~~~~~~~~
/path/to/gnuastro/source/lib/gnuastro/polygon.h:115:1: warning: ‘compareA’ declared ‘static’ but never defined [-Wunused-function]
  115 | compareA(const void *a, const void *b);
      | ^~~~~~~~
/path/to/gnuastro/source/lib/gnuastro/polygon.h:118:1: warning: ‘compareB’ declared ‘static’ but never defined [-Wunused-function]
  118 | compareB(const void *a, const void *b);
      | ^~~~~~~~
/bin/sh ../../libtool  --tag=CC   --mode=link gcc  -Wall -g -O0  -pthread -L\../../lib  -o astcrop main.o ui.o crop.o wcsmode.o onecrop.o ../../bootstrapped/lib/libgnu.la -lgnuastro  -lpthread
libtool: link: gcc -Wall -g -O0 -pthread -o astcrop main.o ui.o crop.o wcsmode.o onecrop.o  -L../../lib ../../bootstrapped/lib/.libs/libgnu.a /usr/local/lib/libgnuastro.so -lgit2 -ltiff -llzma -ljpeg -lwcs -lcfitsio -lcurl -lz -lgsl -lgslcblas -lm -lc -lpthread -pthread -Wl,-rpath -Wl,/usr/local/lib -Wl,-rpath -Wl,/usr/local/lib
/usr/bin/ld: onecrop.o: in function `polygonmask':
/path/to/gnuastro/source/bin/crop/onecrop.c:369: undefined reference to `gal_polygon_vertices_sort_convex'
collect2: error: ld returned 1 exit status
make[2]: *** [Makefile:2001: astcrop] Error 1
make[2]: Target 'all' not remade because of errors.
make[2]: Leaving directory '/dev/shm/source-0.11.155-4d3d/bin/crop'


Didn't you get this error when compiling?

The reason is that you have declared static functions in the `polygon.h' header!!! Static functions are by definition only available within one particular compilation. You can't link to them in other files.

To fix this, remove the declaration of static files from `polygon.h'. Static files must not be in header files (which are loaded into every `.c' file that they are used in).

Mohammad Akhlaghi <makhlaghi>
Group administrator
Thu 02 Apr 2020 12:43:37 PM UTC, comment #36: 

I'll remove whitespaces in the next commit.

Also, there was a warning to include static in header files but not defines? Should I remove these static definitions from the header file. It'll make sense as these are only required in the polygon.c only. What do you think?

Sachin Kumar Singh <sks_15>
Thu 02 Apr 2020 10:13:41 AM UTC, comment #35: 

The functions for concave sort are merged into polygon library and pushed here. Documentation and rebasing will be done before the final push and review. Also, can you please push changes on savannah on GitHub master branch too. It helps to keep track of the number of commits and promotes uniformity:-)

Sachin Kumar Singh <sks_15>
Sun 29 Mar 2020 07:31:21 PM UTC, comment #34: 

I think so! We have all the other necessary components in place now to use it in Crop.

Please add it as a function in the polygon library and update/rebase that branch. Maybe you can call your new function `gal_polygon_vertices_sort'. There is no problem if needs auxiliary functions, they can just be static functions, with the `polygon_' prefix (not `gal_polygon_'). If it is a little long, please group this function and its auxiliary functions in a special group (with 20 empty lines before and after them).

We will also need to rename the old `gal_polygon_ordered_corners' to `gal_polygon_vertices_sort_convex'.

Thanks again for deleting the merged branches and pushing it to your server, it really helps me when viewing the full history on my computer ;-).

Mohammad Akhlaghi <makhlaghi>
Group administrator
Sun 29 Mar 2020 07:18:40 PM UTC, comment #33: 

Is my sorting algorithm for concave polygons ontask-concave-sortready to be used?

Sachin Kumar Singh <sks_15>
Sat 28 Mar 2020 06:00:38 PM UTC, comment #32: 

Previous comment was actually for task #15578, I mistakenly wrote it here! Please reply to it there ;-)

Mohammad Akhlaghi <makhlaghi>
Group administrator
Sat 28 Mar 2020 05:40:47 PM UTC, comment #31: 

I was able to reproduce the problem with the segmentation fault when setting `i' to `size_t'. The problem is in the conditional `i>=0': `size_t' (which is actually `unsigned long') will always be positive, so this test will always pass and it will just read into un-allocated memory until it crashes!

This is an important point when working with unsigned values. To avoid such problems, when looping over unsigned values, always start from `0' and go up. Don't start from the end and come down. You should be able to fix the problem by following this strategy ;-).

Another small tip with the `tests/during-dev.sh' script: before starting its operations, the script will `cd' into `outdir'. So if your testing dataset is in the same directory, you don't need to add its full directory again, just call the file directly ;-). This helps in readability.

Also, since your testers will generally be involved in other things at the moment they are testing, please give some context to what they should expect (what debugging info will be printed), so they don't have to spend time reading the code of your checks to understand what you were doing ;-).

Before checking `gal_polygon_to_counterclockwise', I wanted to see if the lower-level `gal_polygon_iscounterclockwise' works correctly. So I commented your checks in `onecrop.c' temporarily and added these lines to see if it can recognize if the input polygon is sorted or not:


  if(gal_polygon_iscounterclockwise(crp->ipolygon, crp->p->nvertices))
    printf("IS SORTED\n");
  else printf("NOT SORTED\n");


I then added these lines in `tests/during-dev.sh' (under your settings so it overwrites them) to test a simple convex polygon:


# Good convex:
options="--mode=img --polygon=1229,773:2565,1501:1263,2127:823,1304 --polygonnosort"


The output was as expected: it printed `IS SORTED' in the middle of Crop's `stdin' output. The cropped image also looked good in DS9. Then I copied the previous `options', and pasted it under it (to keep the good options for testing later) but let it be overwritten for another test. This time, I just flipped two of the vertices so it becomes a complex polygon:


# Good convex:
options="--mode=img --polygon=1229,773:2565,1501:1263,2127:823,1304 --polygonnosort"

# Bad convex:
options="--mode=img --polygon=1229,773:1263,2127:2565,1501:823,1304 --polygonnosort"


However, it still printed `IS SORTED'! Infact, in the cropped image, you can also see the complex polygon shape (which is good! Showing that your `gal_polygon_isinside' function even works for complex polygons).

But I was expecting it to print `NOT SORTED'.

Can you please see what the problem is? Also, in a similar fashion, add multiple example polygons for various scenarios in `tests/during-dev.sh' to help me test it fast is well ;-).

Mohammad Akhlaghi <makhlaghi>
Group administrator
Sat 21 Mar 2020 10:12:24 PM UTC, comment #30: 

It has been merged with the master branch and pushed to them main repository :-).

After inspecting your commit, I done some modifications to make it blend in better. Then I improved the documentation and finally added a new option in Crop to use these new functions.

Please have a look at the commit messages and the corresponding updates in the code, so we can hopefully have a smoother merge in the future.

I also noticed that you hadn't pulled the updates in the `master' branch before you do the rebase. So a parallel branch was created. This is fine, but generally, when you do a rebase, its best to first do a `git pull' on the master branch to import all recent commits, then rebase it ;-).

With the new `--polygonnosort' option that was just added, Crop can indeed work on Concave polygons, but generally, let's complete task #15578, and also this task (the sorting f concave polygon vertices), then close bug #50120.

Indeed, please mention all those tasks and bugs in your proposal, they are directly related ;-).

Mohammad Akhlaghi <makhlaghi>
Group administrator
Sat 21 Mar 2020 10:52:10 AM UTC, comment #29: 

The requirements are done and the code is rebased and pushed here. Also, as this task fixes bug #50120 and core implementation of task #14639, can I use them in my proposal:-)

Sachin Kumar Singh <sks_15>
Thu 19 Mar 2020 01:31:06 PM UTC, comment #28: 

I just noticed that you haven't actually used your new `gal_polygon_isconvex' function. If you like, you can use it inside of the `polygonmask' of `onecrop.c' and determine the function to use for finding if a given pixel is within the polygon or not.

Since both the `gal_polygon_isinside_convex' and `gal_polygon_isinside' functions have the same types of arguments and return types, you specify a function pointer to set in that conditional and just call that function pointer inside the `POLYGON_MASK' macro ;-).

About the added option part, there is no problem I can do that after I merge your work.

But indeed, I agree, now that you have committed work in Gnuastro's source, you are ready to write a good proposal, so its good to focus these last days on that ;-).

Mohammad Akhlaghi <makhlaghi>
Group administrator
Thu 19 Mar 2020 10:25:01 AM UTC, comment #27: 

I've pushed the code here. I noted that one change in onecrop.c is taking more than 79 columns. If all is okay with other changes I will fix this and rebase and push again. Please review.

I will start working on adding the option soon but I'm a bit preoccupied on my proposal. I've also signed and sent the copyrights form and waiting for a reply.

Sachin Kumar Singh <sks_15>
Wed 18 Mar 2020 09:42:41 PM UTC, comment #26: 

I read the paper about the internal algorithm used in astronomy.net. I understood the algorithm and the required thresholds. I realised how this was related to the non-linear warping task. This is a rough summary:

The astronomy.net using a geometrical hashing algorithm by using predefined indexes from an existing survey and uses a quad asterism. Thereafter they divide the catalogue to be searched on in various healpixel and find the required index after a few passes and if the result passes a bayesian threshold. It used a 4 kd-tree to store and retrieve hashes. Now finally the wcs info of the index is used in the calculation of the positions and any transformations. Also wcs information in the fits files can store any distortions information in them(the key-value pairs)

But now if such a system already exists(and already written in c), what should I do about it for the summer project. Also, the info in all those papers is a little overwhelming at first. Hence please tell me a brief overview of how should I proceed to make them all fit together to make up my project which can be completed within this span of time and is not too complex to implement?



Sachin Kumar Singh <sks_15>
Wed 18 Mar 2020 03:31:04 PM UTC, comment #25: 

Great! I am happy the problem is solved :-).

Here are some points:

0. After a branch is merged into master, please delete it on your remote server also, so it also goes away from my list of branches when I "fetch" your work with `-p'. This is mainly about your old branch that was already merged, but also soon for this one ;-).
0. I noticed that you had committed these two functions over your old `task-concave-sort' branch. I can't go into the details to see if it has all been removed or not. But please make sure that in the rebase we only have the two `gal_polygon_isconvex' and `gal_concavesort_isinside_polygon' functions. We'll get to the sorting problem after this.
0. When doing a `git log -p', I see several lines that end with a white-space (marked in red by Git). It helps if you set this as a rule of thumb: after every commit, run `git log -p -1' and remove such cases. If you use an advanced editor (like GNU Emacs), you can activate a mode to automatically remove such things.

  1. About the naming of the functions, the `gal_polygon_isconvex' is good. But the `gal_polygon_isinside_polygon' function needs to be corrected, since it is a more general implementation of the `gal_polygon_pin' function ("pin" was suppose to mean "point-in"). So my recommendation is to change the name of `gal_polygon_pin' to `gal_polygon_isinside_convex' and your new `gal_polygon_isinside_polygon' can be called `gal_polygon_isinside' (a more general scenario). Also correct all the calls to `gal_polygon_pin' throughout Gnuastro: in the documentation, in the other programs and libraries (a simple Grep will do the job of finding them for you ;-)).
  2. In the book, please consider that we follow the convention of one-sentence-per-line (to help with Git commits). It hasn't yet been fully implemented in the library part of the book yet, but you'll see this if you go up to the source of earlier chapters.


I'll give more comments once its rebased (hopefully there won't be too many more). Until the copyright form is finalized, you can also look into adding the option if you like. I have added a summary of the steps in Comment 7 of Task 15317. You can call it something like `--nopolygonsort'.

Mohammad Akhlaghi <makhlaghi>
Group administrator
Wed 18 Mar 2020 01:08:34 PM UTC, comment #24: 

The code is fixed and tested now. I've also added the documentation for the book. Review it and confirm so that I can rebase the commits.
https://github.com/sachinkumarsingh092/gnuastro/tree/task-is-in-polygon

Sachin Kumar Singh <sks_15>
Wed 18 Mar 2020 11:11:48 AM UTC, comment #23: 

Great! I am happy you found the solution ;-).

Indeed, this is a common problem in many polygon scenarios: you need the first/last point two times! Back when I was writing the polygon.c file, found a very elegant solution to avoid an extra element in the loop. It was introduced in the Second edition of "Computational Geometry in C" book by Joseph O'Rourke (if I am not mistaking!).

Note that this loop is done on every pixel in the image, so its effect on the algorithm's efficiency is significant and we should avoid adding an extra element as much as possible.

You can see the solution in the loops of `gal_polygon_area' and `gal_polygon_pin' (and some others). The last commented paragraph of `gal_polygon_area' gives a summary. I hope this helps ;-).

Mohammad Akhlaghi <makhlaghi>
Group administrator
Wed 18 Mar 2020 07:16:23 AM UTC, comment #22: 

Also, I contacted `assign@gnu.org` again today. I also send a mail to `copyright-clerk@fsf.org` a few days back but haven't got a reply yet.

Sachin Kumar Singh <sks_15>
Wed 18 Mar 2020 07:06:04 AM UTC, comment #21: 

Two new functions are added in `polygon.c' and pushed to `task-is-in-polygon' branch :
https://github.com/sachinkumarsingh092/gnuastro/tree/task-is-in-polygon

Also, I figured out the problem with the extra point and fixed it manually, but need to implement it as a general solution. The problem is that my function need n+1 points(the last edge from v[n-1] to v[0] is also required). We can do this by checking the type of polygon. If it is concave, we need to `realloc(ipolygon, 2*(nvertices+1)*sizeof *ipolygon)` the `ipolygon` to insert 2 more points(the first points) and similarly ordinds. But while doing so, error arise:
` ignoring return value of ‘realloc’, declared with attribute warn_unused_result`.

How to implement my function so that it can use the n+1 points?


Sachin Kumar Singh <sks_15>
Wed 18 Mar 2020 02:04:47 AM UTC, comment #20: 

Thanks for the debug! I think we are getting close.

There is a problem with your test though: the loop in `POLYGON_MASK' is not over the number of polygon vertices, its over the number of pixels in the image: `size'. So the reason you saw zeros after the last polygon vertice was that it was reading blank memory (un-allocated space).

To check the polygon points, I simply put your `printf' statement after this line (in `polygonmask' of `onecrop.c'):


ipolygon[i*2+1] = crp->ipolygon[ordinds[i]*2+1] - fpixel_i[1];


The result was then correct (with all the extra lines printed in the middle of Crop's standard output):


Crop 0.11.33-51ee started on Wed Mar 18 01:52:48 2020
  - Read metadata of 1 dataset.                        0.000659 seconds
-0.180500, 90.761900
155.586300, 245.695800
168.914000, 18.292840
400.481900, 142.406600
566.244500, -0.032670
588.734900, 420.621200
395.484000, 433.948800
561.246600, 544.734900
278.034100, 520.578500
68.956700, 408.959500
  ---- test_cropped.fits              created: 1 input.
Crop finished in:  0.031866 seconds


However, given that the shape you showed perfectly fits the wrong result. I think you should check your function again, probably somewhere in the function you are assuming an extra element in the polygon (that initializes to zero), or are implicitly assuming a zero at the end/start.

I think it would be good to move this part of the task to another branch (with only the changes we are talking about now) and implement the function to find the type of the polygon also. So effectively its just adding two new functions to the polygon library. We can then easily add the option to disable sorting of the polygon vertices, and its ready to be merged :-).

Then we can focus on the other problem of sorting concave polygon vertices.

By the way, if you have no news from the GNU copyright issue, please send them another email and tell them that your work is almost ready to be merged and you are just waiting for them.

Mohammad Akhlaghi <makhlaghi>
Group administrator
Tue 17 Mar 2020 04:05:56 PM UTC, comment #19: 

For the function to find if the polygon is convex or concave, we can simply:

sort all the points;
test all three consecutive points in the list using right-hand cross rule and vector cross product;
  if any of the points fail then the polygon is concave;
  else convex;

Sachin Kumar Singh <sks_15>
Tue 17 Mar 2020 03:23:14 PM UTC, comment #18: 

When I put

printf("[%.3f, %.3f]\n", ipolygon[i*2], ipolygon[i*2+1]);

inside the `POLYGON_MASK' macro in onecrop.c inside its `if' statement, after the first 10 correct vertices(in shifted coordinates), it outputs a lot of (0.000, 0.000) after the last correct point ((68.956700, 408.959500) in this case).When I draw the polygon the bend is there due to this (0, 0) coordinate. I think it's a garbage valve. If I try to prevent it using `if' statement, it shows a segmentation fault. How to prevent it?


Sachin Kumar Singh <sks_15>
Tue 17 Mar 2020 02:00:29 PM UTC, comment #17: 

Probably the best way to debug this problem is to actually use it in the Crop program like I described ;-). Once you do that, you'll see that we don't actually do any polygon clipping (which is in the vector-space).

In this context, we are dealing with pixelated images, not a vector space. Polygon clipping is used in Warp for example when we want to map a warped pixel grid onto an original pixel grid.

The first thing that Crop does is to find the bounding box of the polygon, crops it out and also converts the polygon vertices coordinates to the new cropped image coordinates.

It then simply takes each pixel's coordinates and checks if that pixel is inside the polygon or not (with the function you added).

So if you follow those steps and actually put your function in the polygon.c library and call it through Crop (with the manual disabling of sorting the vertices), you can test it in a more real-world scenario ;-).

Mohammad Akhlaghi <makhlaghi>
Group administrator
Tue 17 Mar 2020 12:45:07 PM UTC, comment #16: 

I've tried many points around the problematic region in my `gal_concavesort_isinside_polygon' but there doesn't seem to be any problem at all. But when reproducing the steps I see that the erroneous region still appears. I even applied a rounding error constant as used in `polygon.c`, but the error still persists. I tried a few other points and then there seemed a pattern.

Firstly, the erroneous region is detached from the main polygon(in ds9 it shows a blank pixel value at the apparent joining). I think it is independent of the main polygon.

If the point is outside(as shown in the figure) the intersection of the 2 points is also considered for the cropping region.

Hence I think the problem arises while performing the actual clipping because, during a convex polygon clipping, no extra region is clipped but during a concave clipping errors might arise(if using Sutherland–Hodgman algorithm).

Which algorithm is used in onecrop.c to clip polygons? How to fix if it used a more general clipping algorithm because `gal_concavesort_isinside_polygon' is tested is not showing any problems?




Sachin Kumar Singh <sks_15>
Mon 16 Mar 2020 02:26:39 AM UTC, comment #15: 

Great! To keep things simple, let's focus on the point-in-polygon part (your `gal_concavesort_isinside_polygon' function). I tried it by copying this function (and its `GAL_CONCAVESORT_TRI_CROSS_PRODUCT' macro) in `lib/polygon.c'.

To be similar to the `gal_polygon_pin', I called it `gal_polygon_pin_general' (to mean that its not just concave), later we should call the old one something like `gal_polygon_pin_convex'. The only correction was that in the end, I had to change the return value to `wn!=0' (we want it to return 1 if the point is inside the polygon and 0 if it isn't).

Then, in Crop's source, I replaced the new function with `gal_polygon_pin' (in the `onecrop_ipolygon_fl' function). In the same `onecrop.c' file, I also disabled sorting manually with this statement around the `gal_polygon_ordered_corners' function.


  if(0)
    gal_polygon_ordered_corners(crp->ipolygon, crp->p->nvertices, ordinds);
  else
    for(i=0;i<crp->p->nvertices;++i) ordinds[i]=i;


Then I called crop on this test.fits file with the following command (you can see the polygon in green on the left window of the attached image).


astcrop test.fits --mode=img --polygon=1318.8195,1043.7619:1474.5863,1198.6958:1487.914,971.29284:1719.4819,1095.4066:1885.2445,952.96733:1907.7349,1373.6212:1714.484,1386.9488:1880.2466,1497.7349:1597.0341,1473.5785:1387.9567,1361.9595


The resulting crop did indeed work nicely compared to the original `gal_polygon_pin'. But if you look in the right window o f the attached image, there is only a small problem in the bottom-left.

If you can fix this, we can finalize the merge of `gal_polygon_pin_general' into the master branch.

It would also be very good if you can add a check function in the polygon.h functions to see if the given polygon is convex or not. The reasons is that if we know this before hand, and the polygon is indeed convex, we can call the simpler/faster `gal_polygon_pin', and if its concave, we can call the more general/complex/slower solution.

The only reason I am suggesting this now is that this is independent of the concave polygon sorting problem and I am really curious about your method in solving it and want to study it a little more in the next few days (hoping that I find more free time ;-)).


Mohammad Akhlaghi <makhlaghi>
Group administrator
Sun 15 Mar 2020 06:35:42 PM UTC, comment #14: 

I've added a new function inside the concavesort.c which uses winding number algorithm to identify points inside a non-convex polygon.
https://github.com/sachinkumarsingh092/gnuastro/tree/task-concave-sort

I've also removed the global variables which really were a nuisance and instead used pointers:-). Some minor bugs were fixed and better names are assigned to the functions. This also fixes task #14639 to a certain extent, only proper sanity checks and formatting for the SAO DS9 format remains. 

Points are now sorted in the anti-clockwise direction.
Regarding the option to avoid sorting and let the program use the same order that the user-provided, I think we can take the input as an edge list and use a flag to indicate that the user is giving such an input. Then we can use a different sort function to sort the edges instead of points(which will be easier to implement using a map/tree data structure).

Sachin Kumar Singh <sks_15>
Sun 15 Mar 2020 02:35:26 AM UTC, comment #13: 

After the last post, and based on the taught that there is possibly no unique concave polygon from a set of points. So I thought it is probably necessary to add an option to avoid sorting and let the program use the same order that the user provided.

I did that manually, but then noticed that even when we do solve the sorting problem, we will have another problem: the function `gal_polygon_pin' (which is used to see if a point is within a polygon) also implicitly assumes a convex polygon. Best said by the function's comment: "If the point is inside the polygon, it will always be to the left of the edge connecting the two vertices when the vertices are traversed in order."

I had a fast search and found this StackOverflow question on the subject with an apparently good (and very well explained!) answer.

Of course, this doesn't mean that we should not sort the vertices! That is indeed important in a generic scenario (where we can't assume that they are sorted). I just want to say that after the sorting of vertices is done, we also need to write a generic point-in-polygon function that can be used for concave polygons too.

Mohammad Akhlaghi <makhlaghi>
Group administrator
Sun 15 Mar 2020 12:50:29 AM UTC, comment #12: 

Thanks Sachin, the code is very well formatted now and it also prints the final output properly :-).

Just one small tip for later: try your best to not global variables. Just pass the values as arguments. Global variables are a treasure trove of bugs that are extremely hard to find!

Ok, back to the code. The first think I did was to draw a random Polygon on paper and give its coordinates. Here they are (sorted in anti-clockwise):


(1.5,1) (2.5,1) (3,2) (4,1) (4.5,4) (1,3) (2,3)


I then shuffled them and put them in the `main' function:


int
main(){
  size_t num=7;
  double out[14];
  double in[14] = {4,1,  1,3,  2.5,1,  3,2,  4.5,4,  1.5,1,  2,3};

  gal_concavesort_sort_and_merge(in, num, out);
  for(int i=0; i<num; i++){
    printf("%.3f, %.3f\n", out[i*2], out[i*2+1]);
  }
  return 0;
}


But the output was like this. It starts in another point and is clockwise, but that isn't the problem.


gcc concave-sort.c && ./a.out
1.000, 3.000
4.500, 4.000
4.000, 1.000
3.000, 2.000
2.500, 1.000
2.000, 3.000
1.500, 1.000


The problem is that around the point (2.5,1) the intended polygon's points were ordered like this: (1.5,1) (2.5,1) (3,2). However, in the output the two points around it are: (3,2) (2.5,1) (3,2): note that the left neighbor is different. A similar case occurred around the point (1.5,1).

The output of the function is also a correct concave polygon, but its not the one I intended! This actually got me thinking about the uniqueness of the concave polygon problem. What do you think?

Mohammad Akhlaghi <makhlaghi>
Group administrator
Sat 14 Mar 2020 09:57:02 PM UTC, comment #11: 

All bugs related to concave sort are fixed. I've tested it with more inputs and it shows good results. I've pushed the code to my `task_concave_sort` branch.
https://github.com/sachinkumarsingh092/gnuastro/tree/task-concave-sort
I've yet not removed the main function so that it can be tested further. Also, the documentation is to be done for the book, which I'll do after the review when the code is ready to push:-)
While compiling the concave-sort.c when there was an error I just commented out that specific portion(which was trivial code anyway) because I tested only this specific file and not others and hence didn't build the whole repository. But I'm sure after building the repo those errors would be no more.

On a side note: I send the copyright form to -email is unavailable- but haven't got a reply yet. Did you get any information regarding this?

Sachin Kumar Singh <sks_15>
Sat 14 Mar 2020 05:40:32 PM UTC, comment #10: 

Thanks a lot Sachin. I just tried it without any change (with the command below), but got all zeros as output:


$ gcc test.c && ./a.out
0.000000 0.000000
0.000000 0.000000
0.000000 0.000000
0.000000 0.000000
0.000000 0.000000
0.000000 0.000000
0.000000 0.000000


I am a little too busy to debug it, but can you please have another look at the coding? It may also be good to add a few debugging steps in the middle for now, so when I run it, you can find the problem.

By the way, it is easier if you just make a new Git project for this single C file, later, when its put in Gnuastro, you can just delete the repository ;-).

Mohammad Akhlaghi <makhlaghi>
Group administrator
Thu 12 Mar 2020 09:35:24 PM UTC, comment #9: 

The program to sort the concave polygons is ready. It is still a crude code to be modified when importing to the library.
I've used a different algorithm which uses a diagonal vector between the rightmost and leftmost point. Then this vector is used to group the points which are above and below it.
Now consider each segment in the output polygon. The first edge from the leftmost point to the first point in A can't intersect any segment (because there is no segment yet). As we proceed in order by x-coordinate through the points in A, from each point, the next segment is going to the right, and all previous segments are to the left. Thus, as we go from the leftmost point, through all the points of A, to the rightmost point, we will have no intersections. This guarantees that no intersection of point or edges will take place. The same can be argued for B in descending order.
Then we sort both the arrays based on their x-coordinate, the point above and on the vector in ascending order and the point below in descending order.
Finally merging these two arrays gives an output array in the desired sorted format.

I've tested this code and it works well. But there is still a minor problem when I try to use the sizes of the temporary arrays as global variables. The problem and its temporary workaround is commented in the code. Please suggest some way to remove this bug. Apart from the bug, I think we need to just beautify the code before exporting it to the gal library.

(file #48589)

Sachin Kumar Singh <sks_15>
Tue 10 Mar 2020 02:10:58 AM UTC, comment #8: 

That was a great review of some of the existing algorithms, I enjoyed reading it ;-). I am sure you'll find/implement the best solutions and we'll test it in many real-world scenarios :-).

About the question in comment #16 (on any particular criterion). I can say no ;-). In the polygon.h library, the functions should be general/abstract, so you don't need to make any assumptions. If caller has any limitations on the points (for example if they are within a data set), its their responsibility before calling this function.

About the GSoC proposal, I agree, its a good time to start writing it and have sufficient time to perfect it without rush. Indeed the two tasks and bug you mentioned are directly related to this task. In fact just today, I was talking with a Gnuastro user about how useful it would be to have task #13565 in Gnuastro (non-linear warpings). In particular, for applying the sky curvature on an image of the sky (which is flat). But there are so many other cases that it can be useful. This is necessary when we want to use Gnuastro to add large images of the sky.

On the same line of thought (to complete Gnuastro in terms of reducing raw observations), an important thing that was also mentioned in the GSoC suggestions is image registration  (finding the warp necessary to align one image with another). This is related to task #14131. These are all also pretty important image processing operations that are also critical in astronomical data reduction. I think focusing on these related tasks (steps done on raw exposures) can be a very good GSoC proposal, and useful for MANY astronomers.

Mohammad Akhlaghi <makhlaghi>
Group administrator
Mon 09 Mar 2020 08:12:32 PM UTC, comment #7: 

Also, I think this is the time around which I should start working on my proposal for GSoC. I think this task and task #14639 and bug #50120 are interconnected and solving this will very much solve the others too and hence should be included for the task for summers.

Also, task #13565 seems interesting of which further discussion will be done on its respective page. This also seems to have a lot of functionalities to be implemented in and seems to be a good task.

Can you suggest some other programs to be implemented for summers and any suggestions for my proposal?

Sachin Kumar Singh <sks_15>
Mon 09 Mar 2020 08:00:45 PM UTC, comment #6: 

After a little bit of research, I found that to uniquely sort the concave polygon we need more than just a set of points as discussed in this answer. The polar sort method using centroid discussed earlier is bound to fail in many cases.

If we want the polygon to be constructed using the outermost set of points, we can use this variation of the gift wrapping algorithm. I think this algorithm will do as it uses lines(or splines?) and angles for its sort.

Also, I've had an idea of my own for the sort using polar sort and distances from a reference point, which will be the point(for anti-clockwise sort) the leftmost(min x-coordinate) and lowermost(min y-coordinate). Now we can sort the points using `atan2`(or just use the polar angle and distance for sorting) and for similar polar angels, the one which is farthest away from it is placed first in the sorted array. This works in O(n) time. This image shows an example. This may not be the best edge case but it is minimal to understand the algorithm. The point O is the reference point. The angle between the reference point and line parallel to x-axis passing through the reference point is the polar angle here.

But first I need to know if there is a criterion for the particular polygon (like it will need to cover particular areas in the file or not or any such case).


Sachin Kumar Singh <sks_15>
Sat 07 Mar 2020 03:09:06 AM UTC, comment #5: 

Indeed, we have to find a robust centroid. Maybe the centroid can be chosen from a set of points within the vertices themselve? Generally, this is a good implementation issue which is easier to answer as you test various implementations. So please try different scenarios and see how successful it is with various sets of points.

About the part of the manual you quoted, it is trying to say why "we will never have a concave polygon" in warping scenarios (where a pixel grid is warped). Initially, I wrote this function to use in Warp, so I didn't need concave polygons at all. But later, when I brought it into Crop too and the library in general, concave polygons also became necessary.

Mohammad Akhlaghi <makhlaghi>
Group administrator
Fri 06 Mar 2020 03:12:14 PM UTC, comment #4: 

What if the centroid lies outside the polygon? Will this method will still be valid? If not, I think we can sort the points using polar coordinates and their distance from a fixed point(which can be the point which is the lowermost and rightmost, for example). But I'm not sure about the feasibility of this method. A bit more discussion is required:-)

Also in the manual, it is stated that the problems with concave polygons *"is because we are always going to be calculating the area of the overlap between
a quadrilateral and the pixel grid or the quadrilateral its self. "* (page:433, polygon.h). Will sorting on this centroid based method have any effect cause we will be calculating centroid on the basis of points and not the area?

Sachin Kumar Singh <sks_15>
Thu 05 Mar 2020 11:16:32 PM UTC, comment #3: 

This seems to be a promising idea! It is an important problem in several astronomical scenarios, not just the Crop program (as mentioned in the first post). I have always been wanting to do this, but never found the time.

In Gnuastro such jobs are done in the polygon.h library. As you see, their inputs and outputs are simple arrays of numbers/coordinates (double *) and integers (size_t). So to start with, you can write and test it completely outside of Gnuastro, once its done and we see it works, we can easily copy it into Gnuastro ;-).


Mohammad Akhlaghi <makhlaghi>
Group administrator
Thu 05 Mar 2020 10:42:17 PM UTC, comment #2: 

I have an idea to sort a concave hull based on their centroid. This answer (https://gamedev.stackexchange.com/a/13232) seems to work. What other things do I need to consider before proceeding? What other problems might arise?

Sachin Kumar Singh <sks_15>
Sun 22 Jan 2017 04:50:51 PM UTC, comment #1: 

Also see, bug #50120 (Print an error when polygon is concave).

Mohammad Akhlaghi <makhlaghi>
Group administrator
Thu 25 Jun 2015 02:44:12 AM UTC, original submission:  

Currently the function to sort the vertices of the polygon (orderedpolygoncorners in lib/polygon.c) only works on convex polygons. Because it was originally designed for ImageWarp where we won't be dealing with concave polygons.

Write a function that is able to sort the polygon edges correctly even if the polygon is concave and use that function for the polygon mode of ImageCrop.

Mohammad Akhlaghi <makhlaghi>
Group administrator

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attached Files
file #48614:  Figure_1.png added by sks_15 (31KiB - image/png)
file #48589:  test.c added by sks_15 (5KiB - text/x-csrc)
file #48570:  geogebra_thumbnail.png added by sks_15 (9KiB - image/png)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by sks_15 (Posted a comment)
  • -email is unavailable- added by makhlaghi (Submitted the item)
  •  

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

     

    Follow 11 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2020-04-09 makhlaghi StatusPostponed Done
        Percent Complete70% 100%
        Open/ClosedOpen Closed
    2020-03-21 makhlaghi Percent Complete0% 70%
        Assigned toNone sks_15
    2020-03-17 sks_15 Attached File- Added Figure_1.png, #48614
    2020-03-17 sks_15 Attached File- Added Screenshot@from@2020-03-16@23-08-03.png, #48612
        Attached File- Added Screenshot@from@2020-03-17@00-07-19.png, #48613
    2020-03-16 makhlaghi Attached File- Added point-in-polygon-comment15.png, #48599
    2020-03-12 sks_15 Attached File- Added test.c, #48589
    2020-03-09 sks_15 Attached File- Added geogebra_thumbnail.png, #48570

    Back to the top

    Powered by Savane 3.13-4448.
    Corresponding source code