Added hexapolar and spiral ray patterns

lib/Petzval/Trace.hs

@@ -8,6 +8,9 @@ module Petzval.Trace
   , HitRecord(..)
   , TraceError(..)
   , raytrace
+  -- * Ray patterns
+  , hexapolarPattern
+  , spiralPattern
   ) where
 
 import Linear
@@ -46,14 +49,14 @@ createRay :: (RealFloat a, Mode a, Epsilon a)
           => Maybe a -- ^ The image plane position. If `Nothing`, the object is at infinity
           -> Pupil a -- ^ The entrance pupil to aim at
           -> a       -- ^ Field angle
-          -> (a,a)   -- ^ Normalized pupil coordinates (in the range \([-1,1]\))
+          -> V2 a    -- ^ Normalized pupil coordinates (in the range \([-1,1]\))
           -> Ray a
-createRay (Just objectPlane) Pupil{position=pz,radius=pr} h (px, py) = 
+createRay (Just objectPlane) Pupil{position=pz,radius=pr} h (V2 px py) = 
   Ray source (normalize $ target ^-^ source)
   where dz = pz - objectPlane 
         source = V3 0 (dz * tan h) objectPlane
         target = V3 (px * pr) (py * pr) pz
-createRay Nothing Pupil{position=pz,radius=pr} h (px, py) = Ray source (normalize $ target ^-^ source)
+createRay Nothing Pupil{position=pz,radius=pr} h (V2 px py) = Ray source (normalize $ target ^-^ source)
   where h' = (pi * (-abs h) / 180) -- field angle in rad
         dy = (V3 0 (cos h') (-sin h')) `project` (V3 0 (py * pr) 0)
         dz =  V3 0 (pz * tan h') (pz * cos h')
@@ -124,3 +127,17 @@ raytrace system ray = trace' 1 ray system
           (fray, rest) <- trace' n2 (nray'&_pos._z -~ element ^. thickness) elements
         
           return (fray & _pos._z +~ element ^. thickness, HitRecord { pos=(nray' ^. _pos), opl} : rest)
+
+-- | Spiral pattern. This is somewhat more irregular than the hexapolar pattern. The argument is the number of points
+spiralPattern :: Floating a => Int -> [V2 a]
+spiralPattern n = map (point . fromIntegral) [0..n-1]
+  where npoints = fromIntegral n - 1
+        point n = let r = sqrt (n / npoints)
+                      theta = 2.3999632297286531 * n
+                  in r *^ V2 (sin theta) (cos theta)
+-- | A hexapolar pattern. This is rather optimally distributed; the argument is the number of rings
+hexapolarPattern :: Floating a => Integer -> [V2 a]
+hexapolarPattern nRings = (V2 0 0) : ( mconcat . map ring $ [0..nRings-1])
+  where ring n = let point t = r *^ V2 (sin t) (cos t)
+                     r = (fromIntegral n) / (fromIntegral nRings)
+                 in map (point . (\x -> x / (fromIntegral) n * 3 * pi) . fromIntegral) [0 .. n-1]