lgrcia / prose

Modular image processing pipelines with Python. Built for Astronomy.
https://prose.readthedocs.io
MIT License
54 stars 11 forks source link

Plate Solving not returning a correct solution: more files #125

Open EAlexJ opened 1 year ago

EAlexJ commented 1 year ago

After getting to try the improved implementation with quite alot of new files from the same Instrument, I noticed that not all the files are being solved correctly.

For you to try to replicate my results I have provided all the files in question here. All of these files already have the correct WCS solution in their header.

There are three categories of issues I found:

First, the file labeled no_solve throws the following error:

Traceback (most recent call last):
  File "/Users/alex/Astronomy/CreateProseSolutions.py", line 38, in <module>
    plate.run(image, show_progress=False)
  File "/Users/alex/Astronomy/.venv/lib/python3.11/site-packages/prose/core/sequence.py", line 111, in run
    self._run(loader=loader)
  File "/Users/alex/Astronomy/.venv/lib/python3.11/site-packages/prose/core/sequence.py", line 134, in _run
    block._run(buffer)
  File "/Users/alex/Astronomy/.venv/lib/python3.11/site-packages/prose/core/block.py", line 73, in _run
    self.run(image)
  File "/Users/alex/Astronomy/.venv/lib/python3.11/site-packages/prose/blocks/catalogs.py", line 208, in run
    new_wcs = twirl.compute_wcs(
              ^^^^^^^^^^^^^^^^^^
  File "/Users/alex/Astronomy/.venv/lib/python3.11/site-packages/twirl/__init__.py", line 57, in compute_wcs
    return fit_wcs_from_points(pixel_coords[i].T, SkyCoord(radecs[j], unit="deg"))
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/alex/Astronomy/.venv/lib/python3.11/site-packages/astropy/wcs/utils.py", line 1144, in fit_wcs_from_points
    fit = least_squares(
          ^^^^^^^^^^^^^^
  File "/Users/alex/Astronomy/.venv/lib/python3.11/site-packages/scipy/optimize/_lsq/least_squares.py", line 818, in least_squares
    raise ValueError("`x0` is infeasible.")
ValueError: `x0` is infeasible.

Second, the files I provided in the wrong_solutions folder are ones where the solution returned by twirl is wrong. It has to be noted, that files originating from the same sighting but in different bands do return valid solutions by the looks of it, only these few are for sure not correct.

Third, the last folder no_solution is an entire sighting where the files can be solved by two different solvers, but twirl does not return a valid solution. They are riddled with quite a few artifacts, but it would be nice if twirl and the prose pipeline would be able to handle these as well.

Here is a snippet of how I ran the solver:

    split_location = file_location.split("/")
    my_pixel_scale = 0.16
    if "H" == split_location[7] or "J" == split_location[7] or "K" == split_location[7]: #get the bandpass out of the filename
        my_pixel_scale = 0.3
    myTelescope = Telescope(pixel_scale=my_pixel_scale)
    image = FITSImage(file_location, telescope=myTelescope)
    plate = Sequence(
        [
            blocks.detection.PointSourceDetection(n=30),
            blocks.catalogs.PlateSolve(debug=True),
        ]
    )
    plate.run(image, show_progress=False)

I am using prose 3.2.3

Please let me know if there are any further questions!

-EAlexJ

lgrcia commented 1 year ago

Thanks a lot for the thorough work @EAlexJ!

First of all I notice that most of the time stars in your images are no well detected due to lots of artifacts on the detector or cosmics. To help with that I would recommend using:

blocks.detection.PointSourceDetection(n=30, minor_length=5, min_area=10)

This will filter out most of the spurious detections.

May you try with that and let me know if that helps? I am happy to move forward and have a look at the cases for which this does not help.

EDIT: Just so you know, the exact sequence I used is:

plate = Sequence(
    [
        blocks.Trim((50, 50)),
        blocks.detection.PointSourceDetection(n=30, minor_length=5, min_area=10),
        blocks.catalogs.PlateSolve(debug=True),
    ]
)

As trimming the overscan was necessary to have cleaner stars detection. With that, for the ones in wrong_solution here are the ones I managed to plate solve

I'll try to check the ones with no solution but my guess is that the initial RA-DEC guess is far from the solution (could it be?)

EAlexJ commented 1 year ago

Thank you for taking your time to look at these!

If you visualize the ones tagged no_solution you will see something like this: Raw_Bad

These are images where an error occurred with the telescope. In total, three solvers where run on them, and for some the solutions are all over the place. The other two solvers are agreeing on the images from 2014 in the g, r and H band, so I assume these are correct.

I just wanted to supply you with those, in case that is something worth pursuing to you, more important is the crash and the ones labeled wrong_solution that you already had some success with.

One thing helpful to your general understanding could be, that these images are each from a single sighting, the different bands are recorded at the same time. This means that there should be at least overlap in the images, as they are basically crops from "one" image in different wavelengths.

Eririf commented 8 months ago

"ValueError: x0 is infeasible" means initial out of bound when using lstsq this may occur when input xy<0 since only crpix has finite bounds