Fix bug when empty array passed to AsphericSurface

[mattderstine]

AxisymmetricOpticalSystem passes an empty array of aspheric terms to AsphericSurface. This fix properly handles that case.

[codecov-commenter]

Codecov Report

Merging #323 (dc44f4b) into main (532a810) will decrease coverage by 0.03%.
The diff coverage is 14.28%.

@@            Coverage Diff             @@
##             main     #323      +/-   ##
==========================================
- Coverage   55.11%   55.08%   -0.04%     
==========================================
  Files          79       79              
  Lines        7267     7273       +6     
==========================================
+ Hits         4005     4006       +1     
- Misses       3262     3267       +5     

Impacted Files	Coverage Δ
src/Optical/OpticalSystem.jl	39.25% <0.00%> (-0.63%)	⬇️
src/Geometry/Primitives/AsphericSurface.jl	59.39% <100.00%> (+0.24%)	⬆️

---

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 532a810...dc44f4b. Read the comment docs.

[BrianGun]

The aspherics argument to AsphericLens and AsphericSurface both look like this:

aspherics::Union{Nothing,Vector{Tuple{Int,T}}}

The intent is that the argument should either be a non-empty vector or nothing. Rather than fixing the problem in AsphericSurface it might be more consistent to fix it in the call to AsphericLens from the constructor of AxisymmetricOpticalSystem, line 378 in OpticalSystem.jl:

elseif surface_type == "Aspheric"
                frontaspherics::Vector{Pair{Int,T}}, backaspherics::Vector{Pair{Int,T}} = [
                    [parse(Int, k) => v for (k, v) in params] for params in [frontparams, backparams]
                ]

                newelement = AsphericLens(
                    material, vertices[i-1], frontradius, frontconic, frontaspherics, backradius, backconic,
                    backaspherics, thickness, semidiameter; lastmaterial, nextmaterial, frontsurfacereflectance,
                    backsurfacereflectance
                )()

Presumably one or both of frontaspherics, backaspherics is an empty vector, when it should be nothing. If we decide to go this way then it would make sense to put an assert in the AsphericSurface constructor that throws an exception when passed an empty vector for the aspheric parameters.

An alternative would be to modify the constructor of AsphericSurface so instead of this:

aspherics::Union{Nothing,Vector{Tuple{Int,T}}}

the argument is this:

aspherics::Vector{Tuple{Int,T}}

along with a change in the documentation saying that zero length vector inputs are acceptable. But that would change the interface of this function so we'd need to track down all the calls to AsphericSurface to make sure they were consistent.

@friggog do you have an opinion?

[friggog]

I agree this change should be made in the caller not the constructor definition as it is now

[mattderstine]

@BrianGun, thanks for the suggestion. I'm working on a fix to AxiSymmetricOpticalSystem. But this does illustrate another problem. When the aspheric terms are all zero, then the intersection can be analytically found and the performance significantly improved. However, when execution gets to an AsphericSurface constructor it's too late for that. An assert could be used to flag this mismatch but I'm virtually certain that this would break all sorts of things (like ZernikeSurface). Or do I misunderstand how this is evaluated?

[BrianGun]

@mattderstine Are you talking about the case when the surface is a conic? It looks like Spherical lenses are properly special cased to faster analytic intersection code, but conics are handled by AcceleratedParametricSurface, which will be much slower than an analytic intersection with a conic.

Looks like the proper place to add this is in the AsphericLens function, which is where the special casing for spherical and aspherical surfaces currently is computed.

We could add conics. Do you want to tackle this? I'd love to do it but I'm up to my eyeballs writing lenslet code.

Ray conic surface intersection doesn't look too difficult, here's one reference among many others: http://skuld.bmsc.washington.edu/people/merritt/graphics/quadrics.html.

[mattderstine]

@BrianGun I've removed the test in AsphericSurface and fixed AxisymmetricOpticalSystem to use nothing instead of an empty array. I now get an error that makes no sense to me. Could you have a look at what happens when you run Examples.draw_zoomlenses()? I don't see why the method call doesn't match the method definition.

I can think about adding conics. My own raytrace code already has the equations for this. I follow the method from Welford's book, Aberrations of Optical Systems, but used Mathematica to find equations that look more like the spherical intersection equations that he derived. No promises on timing. I'd like to finish a solution (more of an alternate approach) for the problem with plotting of 2D concave surfaces and then see if I can properly visualize and analyze the cell phone camera lens from Chipman's book. I'm taking the long route to getting to the polarization code by attempting to implement stuff that others can use rather than just do work arounds along the way.

[BrianGun]

@mattderstine is this the error you are seeing?

julia> using OpticSim
Examples.draw_zoomlenses()
julia> Examples.draw_zoomlenses()
ERROR: ArgumentError: reducing over an empty collection is not allowed
Stacktrace:

if so this error is because the array asphericterms is empty. So the call to AsphericSurface is being passed an empty array of aspheric parameters.

In AsphericLens() the code is making it to the part where an AcceleratedParametricSurface is created but the aspherics argument is an empty array:

   if backaspherics !== nothing
            backaspherics = Tuple{Int,S}.(backaspherics)
        end
        surf = AcceleratedParametricSurface(AsphericSurface(semidiameter + backdecenter_l + S(0.01), radius = backradius, conic = backconic, aspherics = backaspherics), nsamples, interface = opticinterface(S, insidematerial, nextmaterial, backsurfacereflectance, interfacemode))

[mattderstine]

That is the error that I was originally fixing. It happens because the code in AsphericSurface isn't expecting an empty array. The error that I get for my latest update to the pull request is:

ERROR: MethodError: no method matching AsphericLens(::OpticSim.GlassCat.Glass, ::Float64, ::Float64, ::Float64, ::Vector{Pair{Int64, Float64}}, ::Float64, ::Float64, ::Type{Nothing}, ::Float64, ::Float64; lastmaterial=Air, nextmaterial=Air, frontsurfacereflectance=0.0, backsurfacereflectance=0.0)
Closest candidates are:
  AsphericLens(::T, ::S, ::S, ::S, ::Union{Nothing, Array{Pair{Int64, S}, 1}}, ::S, ::S, ::Union{Nothing, Array{Pair{Int64, S}, 1}}, ::S, ::S; lastmaterial, nextmaterial, frontsurfacereflectance, backsurfacereflectance, nsamples, frontdecenter, backdecenter, interfacemode) where {T<:OpticSim.GlassCat.AbstractGlass, S<:Real} at /Users/matt/Development/GitHub/OpticSim.jl/src/Optical/Lenses.jl:38

This seems to match to me. The flagged error is on the second Union{} where the call has a ::Type{Nothing}

[BrianGun]

@mattderstine Not sure why I don't have the latest version of your PR. Not enough git fu on my part.

[mattderstine]

I'm just faking it with git so it might be something that I did. I do see the commits(b7277e3) in this PR.

[BrianGun]

@mattderstine found the error. Your function

turnEmptyIntoNothing(list) = length(list)==0 ? Nothing : list

is returning the type Nothing rather than the value nothing, which is the only instance of Nothing. Then Nothing is passed as an argument to the AsphericLens constructor instead of nothing and the type checker complains.

if you do this

turnEmptyIntoNothing(list) = length(list)==0 ? nothing : list

it should fix this bug.

[mattderstine]

@BrianGun Thanks! I now know much more about nothing. You code does fix the error and it works as expected. I have added an assert to the AsphericSurface constructor and I think all is ready for merge.

I have also changed the AxisymmetricLens constructor so that it calls a ConicLens if the aspheric term vector is empty on both surfaces. This has no advantage now, but when ConicLens is upgraded to use a direct solve of the intersection this should be faster.

[friggog]

Can we add a test for this case? Feel like it should have been caught earlier, otherwise looks good

[BrianGun]

@friggog did you mean an assertion to detect if an empty array is passed to AsphericSurface? @mattderstine added this in the latest commit.

[BrianGun]

I'm approving the PR. Thanks @mattderstine, once again for your contribution!