Skip to content
This repository has been archived by the owner on Oct 23, 2022. It is now read-only.

Galran/adding projected paraxial lenses #258

Merged
merged 58 commits into from
Aug 10, 2021

Conversation

galran
Copy link
Collaborator

@galran galran commented Aug 5, 2021

added the HeadEye module inside the ParaxialAnalysis module to support the creation of Paraxial MLAs based on projection onto a Plane, Cylinder or Sphere.

This is a draft version - comments and exports will be added after initial integration to better support the tools and interfaces we will need.

BrianGun added 30 commits July 18, 2021 13:47
Fixes #242

added first set of functions for computing fraction of beam energy which passes from the paraxial display lens through an aperture.
Fixes #242
adding code to compute beamenergy.
Fixes #242

added function to compute virtual point for paraxial lens. Modified the shape types that can be used for ParaxialLens so they share a common distancefrom plane function.
Fixes #242

Added functions for computing area of spherical triangles and convex polygons. Spherical polygon area computation does not appear to be correct.
Fixes #242

spherical area for spherical triangles and spherical polygons appears to work. Takes approx. 300ns to compute the area of a four sided spherical polygon. Longer than I would have expected. Optimize later if necessary.
Fixes #242

simplified calculation of area of spherical polygon. now runs about 3x faster 76ns for 4 sided polygon whereas before it was 250ns.
Fixes #242

added functions to project pupil vertices onto lens plane, to intersect the projected pupil polygon with the lens aperture polygon, and to compute the beam energy.
…play.

added activebox function which computes which pixels can possibly illuminate a planar polygon located in world space. Used to minimize number of pixels for which have to compute beamenergy.
…tionality needed and has many other features that don't make sense in the paraxial case.
Fixes #242
adding display to lensletassembly.
Fixes #242

moved interface, distancefromplane, normal from all of the surfaces that inherit from PlanarShapes into a generic function defined for PlanarShapes.
added transforms to objects in LesletAssembly so points can be put in the correct coordinate frame.
Fixes #242

wrote test programes for virtual point and projection. Moved beam energy code into lenslet assembly.
Fixes #242

fixed error with missing beamenergy.jl file
Fixes #242

modified planar shapes to return vertices in an SMatrix instead of as SVector{SVector} because LazySet.jl only works with SMatrix and Vector.

Still need to change ConvexPolygon to this form.

started work on testintersection function.
Fixes #242

modified ConvexPolygon to represent polygon points as a matrix. This data structure is more compatible with LazySets which is used to many geometric calculations.
…o represent polygon points. Passed unit tests.
Fixes #242

added * operator for points represented as matrices for better compatibility with LazySets.jl. Unfortunately Julia compiler wasn't smart enough to compile efficient code using high level forms so had to write out the matrix multply loop.
BrianGun and others added 7 commits August 4, 2021 11:07
Fixes #242

fixed Ellipse, Hexagon, Rectangle so they all have zero allocations for vertices3D
Fixes #242

fixed Ellipse, Hexagon, Rectangle so they all have zero allocations for vertices3D
…he T and N in the declaration after encountering some problems in the way ParaxialLens is using the ConvexPolygon
… the creation of Paraxial MLAs based on projection onto a Plane, Cylinder or Sphere
…enses

# Conflicts:
#	src/Geometry/Primitives/NonCSG/ConvexPolygon.jl
@codecov-commenter
Copy link

codecov-commenter commented Aug 5, 2021

Codecov Report

Merging #258 (80a4037) into main (776e788) will decrease coverage by 2.02%.
The diff coverage is 0.00%.

Impacted file tree graph

@@            Coverage Diff             @@
##             main     #258      +/-   ##
==========================================
- Coverage   54.97%   52.95%   -2.03%     
==========================================
  Files          72       77       +5     
  Lines        6708     6965     +257     
==========================================
  Hits         3688     3688              
- Misses       3020     3277     +257     
Impacted Files Coverage Δ
src/Geometry/Primitives/NonCSG/ConvexPolygon.jl 61.90% <0.00%> (ø)
src/Geometry/Primitives/NonCSG/PlanarShape.jl 75.00% <ø> (ø)
src/Optical/Paraxial.jl 69.56% <ø> (ø)
...tical/ParaxialAnalysis/HeadEyeModel/Arrangement.jl 0.00% <0.00%> (ø)
.../Optical/ParaxialAnalysis/HeadEyeModel/Examples.jl 0.00% <0.00%> (ø)
...al/ParaxialAnalysis/HeadEyeModel/HeadEyeClasses.jl 0.00% <0.00%> (ø)
src/Optical/ParaxialAnalysis/HeadEyeModel/Misc.jl 0.00% <0.00%> (ø)
src/Vis/HeadEyeModel.jl 0.00% <0.00%> (ø)
... and 3 more

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 776e788...80a4037. Read the comment docs.

sum = Base.sum(points)
centroid = Vec3(sum * (1.0 / (length(points))))

# Calc full 3x3 covariance matrix, excluding symmetries:
Copy link
Contributor

@BrianGun BrianGun Aug 5, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

there is a function Statistics.cov which I think computes the same thing. Would be simpler to use it, unless you are doing something non-standard. We already use the Statistics package so this wouldn't add a dependency.

Might work more efficiently if points were represented as an SMatrix{3,N}. Otherwise it looks like you'd have to convert from a Vector{SVector} to SMatrix.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

will do.

dir = Vec3(xz*yz - xy*zz, det_y, xy*xz - yz*xx)
else
dir = Vec3(xy*yz - xz*yy,xy*xz - yz*xx,det_z)
end
Copy link
Contributor

@BrianGun BrianGun Aug 5, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is this computation time sensitive? If not then it might be simpler to use the SVD on a matrix of the points and use the singular vectors as the coordinate frame of the best fitting plane. This would only be a few lines of code:

y= x .- mean(x,dims=2)

U,_,_ = svd(y)

then choose which of the U vectors you want to align with the x,y,z axes. Or just return the smallest singular vector as the normal to the plane and compute the other two axes elsewhere.

Built in svd on 3x100 Matrix takes 10usec. If you use MMatrix instead it takes 6us. Or you can use my JacobiSVD (available on request) on MMatrix and it takes 3usec.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

i'm not sure. i had some stability issues in the past with SVG in extreme cases - not sure if because of numerical errors or other reason. will check it out.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@galran if you have stability issues then try my JacobiSVD code. It is supposed to compute small singular vectors to higher accuracy than the more common algorithms. And it is twice as fast as the Julia built in SVD.

c = centroid(poly)

l2w = local2world(poly.local_frame)
len = length(poly.local_points)
len = size(poly.local_points)[2]
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is this right? local_points is a Vector{SVector{2,T}}.

size(poly.local_points)[2]

gives an indexing error on my machine.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

will check.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

poly.local_points is an SMatrix - that is why i changed it from length to size

center(h::Shape) = h._center
points(h::Shape) = h._points
coordinates(h::Shape) = h._coordinates

Copy link
Contributor

@BrianGun BrianGun Aug 5, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is more idiomatic Julia and is simpler

getshapes(::Type{Any};blah blah) = @error "Unknown Type ...."

function get_shapes(::Type{Hexagon}, resolution::Tuple{Int,Int}=(2,2), radius=1.5)::Vector{Shape{2, Float64}}
    cells = Repeat.hexcellsinbox(resolution[1],resolution[2]) 
    hexbasis = Repeat.HexBasis1()
    basic_tile = Repeat.tilevertices(hexbasis) * radius

    res = Vector{Shape{2, Float64}}(undef, 0)
    for c in cells
        center = SVector(hexbasis[c[1], c[2]]) * radius
        points = [(SVector(p...) + center) for p in eachrow(basic_tile)]
        center = center
        push!(res, Shape(center, points, c))
    end
    return res
end

function get_shapes(::Type{Rectangle}, resolution::Tuple{Int,Int}=(2,2), size=1.5)::Vector{Shape{2, Float64}}
    if (typeof(size) == Float64)
        w = h = size
    elseif (typeof(size) == Tuple{Float64, Float64})
        w, h = size
    else
        @error "Unsupported size format [$size] - can be either a Float64 for a square or (Float64, Float64) for a rectangle"
    end

    cells = [[c for c in Iterators.product(-resolution[1]:resolution[1], -resolution[2]:resolution[2])]...]
    # hexbasis = Repeat.HexBasis1()
    # basic_tile = Repeat.tilevertices(hexbasis) * radius
    basic_tile = [w h; w -h; -w -h; -w h]

    res = Vector{Shape{2, Float64}}(undef, 0)
    for c in cells
        center = SVector(Float64(c[1])*w*2.0, Float64(c[2])*h*2.0)
        points = [(SVector(p...) + center) for p in eachrow(basic_tile)]
        center = center
        push!(res, Shape(center, points, c))
    end
    return res
end

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

agreed.



function project(shapes::Vector{Shape{2, T}}, csg::OpticSim.CSGTree{T})::Vector{Shape{3, T}} where {T<:Real}
ray_origin_distance = -10.0
Copy link
Contributor

@BrianGun BrianGun Aug 5, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

should ray_origin_distance be an argument to this function? Or does it make sense for it to always be the same? I haven't read the code carefully enough to figure this out. Just asking, not a criticism.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

well, the assumption is that the original shapes are on the XY plane and the CSG surface to project is in the positive Z half space (above the XY plane), so it is safe to shoot from a negative Z position upwards.
IF the CSG surface is complex and moves into the negative Z it might present a problem.

@@ -21,7 +21,7 @@ The local frame defines the plane (spans by the right and up vectors) with the p
the local_polygon_points are given with respect to the local frame and are 2D points.
NOTE: This class uses static vectors to hold the points which will lead to more efficient performance, but should not be used with polygons with more than 20-30 points.
"""
struct ConvexPolygon{N,T<:Real} <: PlanarShape{T}
struct ConvexPolygon{T<:Real, N} <: PlanarShape{T}
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

was there some reason you changed the order of T and N? I did it this way so it would be similar to the StaticArrays type declarations. Does having N first cause problems with the code?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

reversing this change after our conversation.

Copy link
Contributor

@BrianGun BrianGun left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@galran could you look at the comments on the individual files? github is nagging me to leave this comment for you.

Base automatically changed from BrianGun/issue242 to main August 5, 2021 21:32
Ran Gal added 2 commits August 5, 2021 21:10
…f the system

also removed a couple of commented functions
changed the get_shapes function to rely on the Type instead of a symbol
Added an example function:
OpticSim.ParaxialAnalysis.HeadEye.Examples.example1()
BrianGun
BrianGun previously approved these changes Aug 6, 2021
# Conflicts:
#	src/Geometry/Primitives/NonCSG/ConvexPolygon.jl
#	src/Geometry/Primitives/NonCSG/PlanarShape.jl
#	src/Optical/ParaxialAnalysis/ParaxialAnalysis.jl
BrianGun
BrianGun previously approved these changes Aug 9, 2021
added license header
@BrianGun BrianGun merged commit 82e2cf4 into main Aug 10, 2021
@BrianGun BrianGun deleted the galran/addingProjectedParaxielLenses branch August 10, 2021 02:40
Sign up for free to subscribe to this conversation on GitHub. Already have an account? Sign in.
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants