Add ellipse fitting example

[langestefan]

Draft PR for the ellipse fitting example, will add code.

For this demo I took an image of the cartwheel galaxy captured by the James Webb space telescope:
James Webb Space Telescope images

It's a very general setup and multiple ellipse model definitions, optimization criteria and demo images can be defined.

To efficiently fit the ellipses I need to do some image preprocessing, some filtering and some clustering. For these tasks I used Clustering.jl, Images.jl and DSP.jl. Would it be okay to use these dependencies? I can come up with alternative implementations but it will be more messy.

After fitting the ellipses the final result looks like this:

[odow]

This looks like a super cool (and visual) example.

No problem adding the extra packages. Just ensure they have a compat:

JuMP.jl/docs/Project.toml

Line 42 in 488a0db

[compat]

Also: we need to make sure that the image we use can be redistributed. If possible, find the original source, instead of new scientist?

[langestefan]

Yeah JWST really creates some breathtaking images!

Also: we need to make sure that the image we use can be redistributed. If possible, find the original source, instead of new scientist?

Found the source: https://webbtelescope.org/contents/media/images/2022/039/01G8JXN0K2VBQP112RNSQWTCTH

Unless otherwise stated, all material on the site was produced by NASA and the Space Telescope Science Institute (STScI). Material on this site produced by STScI was created, authored, and/or prepared for NASA under Contract NAS5-03127. Unless otherwise specifically stated, no claim to copyright is being asserted by STScI and material on this site may be freely used as in the public domain in accordance with NASA's contract. However, it is requested that in any subsequent use of this work NASA and STScI be given appropriate acknowledgement.

I'll make sure we add the acknowledgment.

[langestefan]

@odow I tested it locally, it should run.

I have tried to follow the style guide based on other examples but might have missed something.

Maybe good to check if all the MOI calls are correct, and you will probably spot some improvements in JuMP syntax as well.

The final output should look like:

[codecov]

Codecov Report

All modified and coverable lines are covered by tests ✅

Project coverage is 99.58%. Comparing base (00c3b35) to head (8e8e54b).
Report is 2 commits behind head on master.

Additional details and impacted files

@@           Coverage Diff           @@
##           master    #3895   +/-   ##
=======================================
  Coverage   99.58%   99.58%           
=======================================
  Files          43       43           
  Lines        6043     6043           
=======================================
  Hits         6018     6018           
  Misses         25       25           

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

[langestefan]

All comments should be addressed now.

[langestefan]

@odow Thanks for pushing it across the finishline! Seems like it's almost ready, let me know if I can do anything else.

Just curious, but why did the solver change to Clarabel?

[odow]

I changed to Clarabel because SCS wasn't solving the final problem; it was hitting an iteration limit. (As an aside, you must always check the termination condition of a solver before accessing the solution.) But it seems like even Clarabel has trouble.

[langestefan]

You can try changing epsilon from 1e-5. This is the numerical issue with the SPD matrix constraint I was talking about on slack.

[odow]

What happens if we just drop the constraint? Then you might end up with an ellipse with zero area?

[langestefan]

Mathematically that doesn't make sense I think: https://en.m.wikipedia.org/wiki/Matrix_representation_of_conic_sections

I believe I tried this and it just produced trivial forms of the quadratic equation.

I think what we should do is fix the scale to 1. This means we constrain the determinant of Q to be 1.

There is probably a way to do that with MOI. Maybe MOI.logdetconetriangle or similar? I'll check tomorrow if you haven't already (I'm in EUW timezone ;) )

[langestefan]

I just checked the results with the SPD constraint removed.

Using SCS:

Using Clarabel:

This is quite puzzling to me. I am not suprised SCS now fails, but that Clarabel succeeds I don't understand.

[odow]

This is quite puzzling to me. I am not suprised SCS now fails

It's not that SCS now fails. It didn't work before (it hit an iteration limit). It was just luck that the primal solution that SCS terminated with looked okay.

but that Clarabel succeeds I don't understand.

Different solvers different algorithms 😄

[odow]

Thanks! That took a while to get right, but we got there in the end 😄

[langestefan]

Yay! Thanks for finishing it. I learned a lot again.