ModuleMirror
Provides classes for different mirror surfaces, which are types of optics.
Think of these as the z-coordinates on top of the x-y-grid provided by the ART.ModuleSupport.Support.
Also provides the function ReflectionMirrorRayList that returns the rays reflected on the mirror. Rays that do not hit the support are not propagated any further.
Created in 2019
@author: Anthony Guillaume and Stefan Haessler and Andre Kalouguine and Charles Bourassin-Bouchet
1""" 2Provides classes for different mirror surfaces, which are types of optics. 3 4Think of these as the z-coordinates on top of the x-y-grid provided by the ART.ModuleSupport.Support. 5 6Also provides the function *ReflectionMirrorRayList* that returns the rays reflected on 7the mirror. Rays that do not hit the support are not propagated any further. 8 9 10 11 12 13Created in 2019 14 15@author: Anthony Guillaume and Stefan Haessler and Andre Kalouguine and Charles Bourassin-Bouchet 16""" 17# %% Modules 18import ARTcore.ModuleGeometry as mgeo 19from ARTcore.ModuleGeometry import Point, Vector, Origin 20import ARTcore.ModuleSupport as msup 21import ARTcore.ModuleOpticalRay as mray 22import ARTcore.ModuleOpticalElement as moe 23from ARTcore.DepGraphDefinitions import OAP_calculator, EllipsoidalMirrorCalculator 24import ARTcore.ModuleSurface as msurf 25 26import numpy as np 27from enum import Enum 28import math 29from abc import ABC, abstractmethod 30from copy import copy 31import logging 32import time 33 34logger = logging.getLogger(__name__) 35 36# %% Generic mirror class definition 37 38class Curvature(Enum): 39 CONVEX = -1 40 FLAT = 0 41 CONCAVE = 1 42 43class Mirror(moe.OpticalElement): 44 """Abstract base class for mirrors.""" 45 vectorised = False 46 def __init__(self): 47 self.Surface = msurf.IdealSurface() 48 super().__init__() 49 50 @abstractmethod 51 def _get_intersection(self, Ray): 52 """ 53 This method should return the intersection point between the ray and the mirror surface. 54 The calculation should be done in the reference frame of the mirror. 55 It essentially defines the mirror surface (but not its normal). 56 """ 57 pass 58 59 @abstractmethod 60 def get_local_normal(self, Point: np.ndarray): 61 """ 62 This method should return the normal unit vector in point 'Point' on the mirror surface. 63 The normal is in the reference frame of the mirror. 64 """ 65 pass 66 67 def propagate_raylist(self, RayList, alignment=False): 68 """ 69 Propagate a list of rays to the mirror and return the reflected rays. 70 71 Parameters 72 ---------- 73 RayList : a RayList (a class defined in ARTcore.ModuleOpticalRay) 74 Rays to propagate through the mirror. 75 76 alignment : bool 77 If True, the rays are not blocked by the support. Not implemented yet 78 79 Returns 80 ------- 81 ReflectedRayList : a RayList 82 The rays that were reflected by the mirror. 83 """ 84 #local_rays = [ray.to_basis(*self.basis) for ray in RayList] 85 local_rays = RayList.to_basis(*self.basis) 86 # So local_rays is a list of rays in the mirror reference frame 87 if self.vectorised and len(local_rays) > 10: 88 # If the class implements the vectorised version of the intersection and normal calculation 89 intersection_points, OK = self._get_intersections(local_rays) 90 # NB: the _get_intersections method should return two things: 91 # - the intersection points in the mirror reference frame with nan's for rays that do not intersect 92 # - a boolean array indicating which rays intersections make sense 93 local_normals = mgeo.VectorArray(np.zeros_like(intersection_points)) 94 local_normals[OK] = self.get_local_normals(intersection_points[OK]) 95 local_reflected_rays = [ 96 self.Surface.reflect_ray(local_rays[i], intersection_points[i], local_normals[i]) 97 for i in range(len(local_rays)) if OK[i]] 98 else: 99 intersection_points, OK = zip(*[self._get_intersection(r) for r in local_rays]) 100 # OK is a boolean list indicating if the ray intersects the mirror (for instance is in the support) 101 local_reflected_rays = [ 102 self.Surface.reflect_ray(local_rays[i], intersection_points[i], self.get_local_normal(intersection_points[i])) 103 for i in range(len(local_rays)) if OK[i]] 104 reflected_rays = mray.RayList.from_list(local_reflected_rays).from_basis(*self.basis) 105 if alignment and len(reflected_rays) == 0: 106 logger.error("Unexpected error occurred: no rays were reflected during alignment procedure.") 107 logger.error("This should not happen.") 108 logger.error(f"Mirror on which error ocurred: {self}") 109 logger.error(f"Alignment ray: {RayList[0]}") 110 logger.error(f"Alignment ray in mirror reference frame: {local_rays[0]}") 111 logger.error(f"Intersection point: {self._get_intersection(local_rays[0])[0]}") 112 logger.error(f"Support: {self.support}") 113 elif len(reflected_rays) == 0: 114 logger.warning("No rays were reflected by the mirror.") 115 logger.debug(f"Mirror: {self}") 116 logger.debug(f"First ray: {RayList[0]}") 117 logger.debug(f"First ray in mirror reference frame: {RayList[0].to_basis(self.r0, self.r, self.q)}") 118 logger.debug(f"First ray intersection point: {self._get_intersection(RayList[0].to_basis(self.r0, self.r, self.q))}") 119 return reflected_rays 120 121 122def _IntersectionRayMirror(PointMirror, ListPointIntersectionMirror): 123 """When the ray has pierced the mirror twice, select the first point, otherwise just keep that one.""" 124 if len(ListPointIntersectionMirror) == 2: 125 return mgeo.ClosestPoint( 126 PointMirror, 127 ListPointIntersectionMirror[0], 128 ListPointIntersectionMirror[1], 129 ) 130 elif len(ListPointIntersectionMirror) == 1: 131 return ListPointIntersectionMirror[0] 132 else: 133 return None 134 135 136# %% Plane mirror definitions 137class MirrorPlane(Mirror): 138 """ 139 A plane mirror. 140 141 Attributes 142 ---------- 143 support : ART.ModuleSupport.Support 144 Mirror support 145 type : str = 'Plane Mirror' 146 Human readable mirror type 147 Methods 148 ------- 149 MirrorPlane.get_normal(Point) 150 151 MirrorPlane.get_centre() 152 153 MirrorPlane.get_grid3D(NbPoints) 154 """ 155 vectorised = True 156 def __init__(self, Support, **kwargs): 157 """ 158 Create a plane mirror on a given Support. 159 160 Parameters 161 ---------- 162 Support : ART.ModuleSupport.Support 163 164 """ 165 super().__init__() 166 self.support = Support 167 self.type = "Plane Mirror" 168 self.curvature = Curvature.FLAT 169 170 if "Surface" in kwargs: 171 self.Surface = kwargs["Surface"] 172 173 self.r0 = Point([0.0, 0.0, 0.0]) 174 175 self.centre_ref = Point([0.0, 0.0, 0.0]) 176 self.support_normal_ref = Vector([0, 0, 1.0]) 177 self.majoraxis_ref = Vector([1.0, 0, 0]) 178 179 self.add_global_points("centre") 180 self.add_global_vectors("support_normal", "majoraxis") 181 182 def _get_intersection(self, Ray): 183 """Return the intersection point between Ray and the xy-plane.""" 184 t = -Ray.point[2] / Ray.vector[2] 185 PointIntersection = Ray.point + Ray.vector * t 186 return PointIntersection, t > 0 and PointIntersection in self.support 187 188 def get_local_normal(self, Point: np.ndarray): 189 """Return normal unit vector in point 'Point' on the plane mirror.""" 190 return Vector([0, 0, 1]) 191 192 def _get_intersections(self, RayList): 193 """ 194 Vectorised version of the intersection calculation. 195 """ 196 vectors = np.array([ray.vector for ray in RayList]) 197 points = np.array([ray.point for ray in RayList]) 198 199 t = -points[:,2] / vectors[:,2] 200 Points = points + vectors * t[:,np.newaxis] 201 OK = (t > 0) & np.array([p in self.support for p in Points]) 202 return mgeo.PointArray(Points), OK 203 204 def get_local_normals(self, Points): 205 return mgeo.VectorArray(np.zeros_like(Points) + [0, 0, 1]) 206 207 def _zfunc(self, PointArray): 208 x = PointArray[:,0] 209 y = PointArray[:,1] 210 return np.zeros_like(x) 211 212# %% Spherical mirror definitions 213class MirrorSpherical(Mirror): 214 """ 215 Spherical mirror surface with eqn. $x^2 + y^2 + z^2 = R^2$, where $R$ is the radius. 216 217  218 219 Attributes 220 ---------- 221 radius : float 222 Radius of curvature. A postive value for concave mirror, a negative value for a convex mirror. 223 224 support : ART.ModuleSupport.Support 225 226 type : str SphericalCC Mirror' or SphericalCX Mirror' 227 228 Methods 229 ------- 230 MirrorSpherical.get_normal(Point) 231 232 MirrorSpherical.get_centre() 233 234 MirrorSpherical.get_grid3D(NbPoints) 235 236 """ 237 def __init__(self, Support, Radius, **kwargs): 238 """ 239 Construct a spherical mirror. 240 241 Parameters 242 ---------- 243 Radius : float 244 The radius of curvature in mm. A postive value for concave mirror, a negative value for a convex mirror. 245 246 Support : ART.ModuleSupport.Support 247 248 """ 249 super().__init__() 250 if Radius < 0: 251 self.type = "SphericalCX Mirror" 252 self.curvature = Curvature.CONVEX 253 self.radius = -Radius 254 elif Radius > 0: 255 self.type = "SphericalCC Mirror" 256 self.curvature = Curvature.CONCAVE 257 self.radius = Radius 258 else: 259 raise ValueError("Radius of curvature must be non-zero") 260 self.support = Support 261 262 if "Surface" in kwargs: 263 self.Surface = kwargs["Surface"] 264 265 self.r0 = Point([0.0, 0, -self.radius]) 266 267 self.centre_ref = Point([0.0, 0, -self.radius]) 268 self.sphere_center_ref = Point([0.0, 0, 0]) 269 self.focus_ref = Point([0.0, 0, -self.radius/2]) 270 271 self.support_normal_ref = Vector([0, 0, 1.0]) 272 self.majoraxis_ref = Vector([1.0, 0, 0]) 273 self.towards_focus_ref = (self.focus_ref - self.centre_ref).normalized() 274 275 self.add_global_points("sphere_center", "focus", "centre") 276 self.add_global_vectors("towards_focus", "majoraxis", "support_normal") 277 278 def _get_intersection(self, Ray): 279 """Return the intersection point between the ray and the sphere.""" 280 a = np.dot(Ray.vector, Ray.vector) 281 b = 2 * np.dot(Ray.vector, Ray.point) 282 c = np.dot(Ray.point, Ray.point) - self.radius**2 283 284 Solution = mgeo.SolverQuadratic(a, b, c) 285 Solution = mgeo.KeepPositiveSolution(Solution) 286 287 ListPointIntersection = [] 288 for t in Solution: 289 Intersect = Ray.vector * t + Ray.point 290 if Intersect[2] < 0 and Intersect in self.support: 291 ListPointIntersection.append(Intersect) 292 293 return _IntersectionRayMirror(Ray.point, ListPointIntersection) 294 295 def get_local_normal(self, Point): 296 """Return the normal unit vector on the spherical surface at point Point.""" 297 return -Vector(Point).normalized() 298 299 def _zfunc(self, PointArray): 300 x = PointArray[:,0] 301 y = PointArray[:,1] 302 return -np.sqrt(self.radius**2 - x**2 - y**2) 303 304# %% Parabolic mirror definitions 305class MirrorParabolic(Mirror): 306 r""" 307 A paraboloid with vertex at the origin $O=[0,0,0]$ and symmetry axis z: 308 $z = \frac{1}{4f}[x^2 + y^2]$ where $f$ is the focal lenght of the *mother* 309 parabola (i.e. measured from its center at $O$ to the focal point $F$). 310 311 The center of the support is shifted along the x-direction by the off-axis distance $x_c$. 312 This leads to an *effective focal length* $f_\mathrm{eff}$, measured from the shifted center 313 of the support $P$ to the focal point $F$. 314 It is related to the mother focal length by $f = f_\\mathrm{eff} \cos^2(\alpha/2) $, 315 or equivalently $ p = 2f = f_\mathrm{eff} (1+\\cos\\alpha)$, where $\\alpha$ 316 is the off-axis angle, and $p = 2f$ is called the semi latus rectum. 317 318 Another useful relationship is that between the off-axis distance and the resulting 319 off-axis angle: $x_c = 2 f \tan(\alpha/2)$. 320 321 322  323 324 Attributes 325 ---------- 326 offaxisangle : float 327 Off-axis angle of the parabola. Modifying this also updates p, keeping feff constant. 328 Attention: The off-axis angle must be *given in degrees*, but is stored and will be *returned in radian* ! 329 330 feff : float 331 Effective focal length of the parabola in mm. Modifying this also updates p, keeping offaxisangle constant. 332 333 p : float 334 Semi latus rectum of the parabola in mm. Modifying this also updates feff, keeping offaxisangle constant. 335 336 support : ART.ModuleSupport.Support 337 338 type : str 'Parabolic Mirror'. 339 340 Methods 341 ------- 342 MirrorParabolic.get_normal(Point) 343 344 MirrorParabolic.get_centre() 345 346 MirrorParabolic.get_grid3D(NbPoints) 347 348 """ 349 #vectorised = True # Unitl I fix the closest solution thing 350 def __init__(self, Support, 351 FocalEffective: float=None, 352 OffAxisAngle: float = None, 353 FocalParent: float = None, 354 RadiusParent: float = None, 355 OffAxisDistance: float = None, 356 MoreThan90: bool = None, 357 **kwargs): 358 """ 359 Initialise a Parabolic mirror. 360 361 Parameters 362 ---------- 363 FocalEffective : float 364 Effective focal length of the parabola in mm. 365 366 OffAxisAngle : float 367 Off-axis angle *in degrees* of the parabola. 368 369 Support : ART.ModuleSupport.Support 370 371 """ 372 super().__init__() 373 self.curvature = Curvature.CONCAVE 374 self.support = Support 375 self.type = "Parabolic Mirror" 376 377 if "Surface" in kwargs: 378 self.Surface = kwargs["Surface"] 379 380 parameter_calculator = OAP_calculator() 381 values,steps = parameter_calculator.calculate_values(verify_consistency=True, 382 fs=FocalEffective, 383 theta=OffAxisAngle, 384 fp=FocalParent, 385 Rc=RadiusParent, 386 OAD=OffAxisDistance, 387 more_than_90=MoreThan90) 388 self._offaxisangle = values["theta"] 389 self._feff = values["fs"] 390 self._p = values["p"] 391 self._offaxisdistance = values["OAD"] 392 self._fparent = values["fp"] 393 self._rparent = values["Rc"] 394 395 self.r0 = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 396 397 self.centre_ref = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 398 self.support_normal_ref = Vector([0, 0, 1.0]) 399 self.majoraxis_ref = Vector([1.0, 0, 0]) 400 401 self.focus_ref = Point([0.0, 0, self._fparent]) 402 self.towards_focusing_ref = (self.focus_ref-self.centre_ref).normalized() 403 self.towards_collimated_ref = Vector([0, 0, 1.0]) 404 405 self.add_global_points("focus", "centre") 406 self.add_global_vectors("towards_focusing", "towards_collimated", "support_normal", "majoraxis") 407 408 def _get_intersection(self, Ray): 409 """Return the intersection point between the ray and the parabola.""" 410 ux, uy, uz = Ray.vector 411 xA, yA, zA = Ray.point 412 413 da = ux**2 + uy**2 414 db = 2 * (ux * xA + uy * yA) - 2 * self._p * uz 415 dc = xA**2 + yA**2 - 2 * self._p * zA 416 417 Solution = mgeo.SolverQuadratic(da, db, dc) 418 Solution = mgeo.KeepPositiveSolution(Solution) 419 IntersectionPoint = [Ray.point + Ray.vector * t for t in Solution if t > 0] 420 distances = [mgeo.Vector(i - self.r0).norm for i in IntersectionPoint] 421 IntersectionPoint = IntersectionPoint[distances.index(min(distances))] 422 return IntersectionPoint, IntersectionPoint-self.r0 in self.support 423 424 def _get_intersections(self, RayList): 425 """ 426 Vectorised version of the intersection calculation. 427 """ 428 vectors = np.array([ray.vector for ray in RayList]) 429 points = np.array([ray.point for ray in RayList]) 430 Solutions = np.zeros(len(RayList), dtype=float) 431 432 da = np.sum(vectors[:,0:2]**2, axis=1) 433 db = 2 * (np.sum(vectors[:,0:2] * points[:,0:2], axis=1) - self._p * vectors[:,2]) 434 dc = np.sum(points[:,0:2]**2, axis=1) - 2 * self._p * points[:,2] 435 436 linear = da == 0 437 Solutions[linear] = -dc[linear] / db[linear] 438 439 centered = (dc==0) & (da!=0) # If equation is of form ax**2 + bx = 0 440 Solutions[centered] = np.maximum(-db[centered] / da[centered], 0) 441 442 unstable = (4*da*dc) < db**2 * 1e-10 # Cases leading to catastrophic cancellation 443 444 Deltas = db**2 - 4 * da * dc 445 OK = Deltas >= 0 446 SolutionsPlus = (-db + np.sqrt(Deltas)) / 2 / da 447 SolutionsMinus = (-db - np.sqrt(Deltas)) / 2 / da 448 SolutionsPlus[unstable] = dc[unstable] / (da[unstable] * SolutionsMinus[unstable]) 449 SolutionsPlus = np.where(SolutionsPlus >= 0, SolutionsPlus, np.inf) 450 SolutionsMinus = np.where(SolutionsMinus >= 0, SolutionsMinus, np.inf) 451 Solutions = np.minimum(SolutionsPlus, SolutionsMinus) 452 Solutions = np.maximum(Solutions, 0) 453 Points = points + vectors * Solutions[:,np.newaxis] 454 OK = OK & (Solutions > 0) & np.array([p-self.r0 in self.support for p in Points]) 455 return mgeo.PointArray(Points), OK 456 457 def get_local_normal(self, Point): 458 """Return the normal unit vector on the paraboloid surface at point Point.""" 459 Gradient = Vector(np.zeros(3)) 460 Gradient[0] = -Point[0] 461 Gradient[1] = -Point[1] 462 Gradient[2] = self._p 463 return Gradient.normalized() 464 465 def get_local_normals(self, Points): 466 """ 467 Vectorised version of the normal calculation. 468 """ 469 Gradients = mgeo.VectorArray(np.zeros_like(Points)) 470 Gradients[:,0] = -Points[:,0] 471 Gradients[:,1] = -Points[:,1] 472 Gradients[:,2] = self._p 473 return Gradients.normalized() 474 475 def _zfunc(self, PointArray): 476 x = PointArray[:,0] 477 y = PointArray[:,1] 478 return (x**2 + y**2) / 2 / self._p 479 480 481# %% Toroidal mirror definitions 482class MirrorToroidal(Mirror): 483 r""" 484 Toroidal mirror surface with eqn.$(\sqrt{x^2+z^2}-R)^2 + y^2 = r^2$ where $R$ and $r$ the major and minor radii. 485 486  487 488 Attributes 489 ---------- 490 majorradius : float 491 Major radius of the toroid in mm. 492 493 minorradius : float 494 Minor radius of the toroid in mm. 495 496 support : ART.ModuleSupport.Support 497 498 type : str 'Toroidal Mirror'. 499 500 Methods 501 ------- 502 MirrorToroidal.get_normal(Point) 503 504 MirrorToroidal.get_centre() 505 506 MirrorToroidal.get_grid3D(NbPoints) 507 508 """ 509 510 def __init__(self, Support, MajorRadius, MinorRadius, **kwargs): 511 """ 512 Construct a toroidal mirror. 513 514 Parameters 515 ---------- 516 MajorRadius : float 517 Major radius of the toroid in mm. 518 519 MinorRadius : float 520 Minor radius of the toroid in mm. 521 522 Support : ART.ModuleSupport.Support 523 524 """ 525 super().__init__() 526 self.support = Support 527 self.type = "Toroidal Mirror" 528 529 if "Surface" in kwargs: 530 self.Surface = kwargs["Surface"] 531 532 self.curvature = Curvature.CONCAVE 533 534 self.majorradius = MajorRadius 535 self.minorradius = MinorRadius 536 537 self.r0 = Point([0.0, 0.0, -MajorRadius - MinorRadius]) 538 539 self.centre_ref = Point([0.0, 0.0, -MajorRadius - MinorRadius]) 540 self.support_normal_ref = Vector([0, 0, 1.0]) 541 self.majoraxis_ref = Vector([1.0, 0, 0]) 542 543 self.add_global_vectors("support_normal", "majoraxis") 544 self.add_global_points("centre") 545 546 def _get_intersection(self, Ray): 547 """Return the intersection point between Ray and the toroidal mirror surface. The intersection points are given in the reference frame of the mirror""" 548 ux = Ray.vector[0] 549 uz = Ray.vector[2] 550 xA = Ray.point[0] 551 zA = Ray.point[2] 552 553 G = 4.0 * self.majorradius**2 * (ux**2 + uz**2) 554 H = 8.0 * self.majorradius**2 * (ux * xA + uz * zA) 555 I = 4.0 * self.majorradius**2 * (xA**2 + zA**2) 556 J = np.dot(Ray.vector, Ray.vector) 557 K = 2.0 * np.dot(Ray.vector, Ray.point) 558 L = ( 559 np.dot(Ray.point, Ray.point) 560 + self.majorradius**2 561 - self.minorradius**2 562 ) 563 564 a = J**2 565 b = 2 * J * K 566 c = 2 * J * L + K**2 - G 567 d = 2 * K * L - H 568 e = L**2 - I 569 570 Solution = mgeo.SolverQuartic(a, b, c, d, e) 571 if len(Solution) == 0: 572 return None, False 573 Solution = [t for t in Solution if t > 0] 574 IntersectionPoint = [Ray.point + Ray.vector *i for i in Solution] 575 distances = [mgeo.Vector(i - self.r0).norm for i in IntersectionPoint] 576 IntersectionPoint = IntersectionPoint[distances.index(min(distances))] 577 OK = IntersectionPoint - self.r0 in self.support and np.abs(IntersectionPoint[2]-self.r0[2]) < self.majorradius 578 return IntersectionPoint, OK 579 580 def get_local_normal(self, Point): 581 """Return the normal unit vector on the toroidal surface at point Point in the reference frame of the mirror""" 582 x = Point[0] 583 y = Point[1] 584 z = Point[2] 585 A = self.majorradius**2 - self.minorradius**2 586 587 Gradient = Vector(np.zeros(3)) 588 Gradient[0] = ( 589 4 * (x**3 + x * y**2 + x * z**2 + x * A) 590 - 8 * x * self.majorradius**2 591 ) 592 Gradient[1] = 4 * (y**3 + y * x**2 + y * z**2 + y * A) 593 Gradient[2] = ( 594 4 * (z**3 + z * x**2 + z * y**2 + z * A) 595 - 8 * z * self.majorradius**2 596 ) 597 598 return -Gradient.normalized() 599 600 def _zfunc(self, PointArray): 601 x = PointArray[:,0] 602 y = PointArray[:,1] 603 return -np.sqrt( 604 (np.sqrt(self.minorradius**2 - y**2) + self.majorradius)**2 - x**2 605 ) 606 607def ReturnOptimalToroidalRadii( 608 Focal: float, AngleIncidence: float 609) -> (float, float): 610 """ 611 Get optimal parameters for a toroidal mirror. 612 613 Useful helper function to get the optimal major and minor radii for a toroidal mirror to achieve a 614 focal length 'Focal' with an angle of incidence 'AngleIncidence' and with vanishing astigmatism. 615 616 Parameters 617 ---------- 618 Focal : float 619 Focal length in mm. 620 621 AngleIncidence : int 622 Angle of incidence in degrees. 623 624 Returns 625 ------- 626 OptimalMajorRadius, OptimalMinorRadius : float, float. 627 """ 628 AngleIncidenceRadian = AngleIncidence * np.pi / 180 629 OptimalMajorRadius = ( 630 2 631 * Focal 632 * (1 / np.cos(AngleIncidenceRadian) - np.cos(AngleIncidenceRadian)) 633 ) 634 OptimalMinorRadius = 2 * Focal * np.cos(AngleIncidenceRadian) 635 return OptimalMajorRadius, OptimalMinorRadius 636 637 638# %% Ellipsoidal mirror definitions 639class MirrorEllipsoidal(Mirror): 640 """ 641 Ellipdoidal mirror surface with eqn. $(x/a)^2 + (y/b)^2 + (z/b)^2 = 1$, where $a$ and $b$ are semi major and semi minor axes. 642 643  644 645 Attributes 646 ---------- 647 a : float 648 Semi major axis of the ellipsoid in mm. 649 650 b : float 651 Semi minor axis of the ellipsoid in mm. 652 653 support : ART.ModuleSupport.Support 654 655 type : str 'Ellipsoidal Mirror'. 656 657 Methods 658 ------- 659 MirrorEllipsoidal.get_normal(Point) 660 661 MirrorEllipsoidal.get_centre() 662 663 MirrorEllipsoidal.get_grid3D(NbPoints) 664 665 """ 666 667 def __init__( 668 self, 669 Support, 670 SemiMajorAxis=None, 671 SemiMinorAxis=None, 672 OffAxisAngle=None, 673 f_object=None, 674 f_image=None, 675 IncidenceAngle=None, 676 Magnification=None, 677 DistanceObjectImage=None, 678 Eccentricity=None, 679 **kwargs, 680 ): 681 """ 682 Generate an ellipsoidal mirror with given parameters. 683 684 The angles are given in degrees but converted to radians for internal calculations. 685 You can specify the semi-major and semi-minor axes, the off-axis angle, the object and image focal distances or some other parameters. 686 The constructor uses a magic calculator to determine the missing parameters. 687 688 Parameters 689 ---------- 690 Support : TYPE 691 ART.ModuleSupport.Support. 692 SemiMajorAxis : float (optional) 693 Semi major axis of the ellipsoid in mm.. 694 SemiMinorAxis : float (optional) 695 Semi minor axis of the ellipsoid in mm.. 696 OffAxisAngle : float (optional) 697 Off-axis angle of the mirror in mm. Defined as the angle at the centre of the mirror between the two foci.. 698 f_object : float (optional) 699 Object focal distance in mm. 700 f_image : float (optional) 701 Image focal distance in mm. 702 IncidenceAngle : float (optional) 703 Angle of incidence in degrees. 704 Magnification : float (optional) 705 Magnification. 706 DistanceObjectImage : float (optional) 707 Distance between object and image in mm. 708 Eccentricity : float (optional) 709 Eccentricity of the ellipsoid. 710 """ 711 super().__init__() 712 self.type = "Ellipsoidal Mirror" 713 self.support = Support 714 715 if "Surface" in kwargs: 716 self.Surface = kwargs["Surface"] 717 718 parameter_calculator = EllipsoidalMirrorCalculator() 719 values, steps = parameter_calculator.calculate_values(verify_consistency=True, 720 f1=f_object, 721 f2=f_image, 722 a=SemiMajorAxis, 723 b=SemiMinorAxis, 724 offset_angle=OffAxisAngle, 725 incidence_angle=IncidenceAngle, 726 m=Magnification, 727 l=DistanceObjectImage, 728 e=Eccentricity) 729 self.a = values["a"] 730 self.b = values["b"] 731 self._offaxisangle = np.deg2rad(values["offset_angle"]) 732 self._f_object = values["f1"] 733 self._f_image = values["f2"] 734 735 self.centre_ref = self._get_centre_ref() 736 737 self.r0 = self.centre_ref 738 739 self.support_normal_ref = Vector([0, 0, 1.0]) 740 self.majoraxis_ref = Vector([1.0, 0, 0]) 741 742 self.f1_ref = Point([-self.a, 0,0]) 743 self.f2_ref = Point([ self.a, 0,0]) 744 self.towards_image_ref = (self.f2_ref - self.centre_ref).normalized() 745 self.towards_object_ref = (self.f1_ref - self.centre_ref).normalized() 746 self.centre_normal_ref = self.get_local_normal(self.centre_ref) 747 748 self.add_global_points("f1", "f2", "centre") 749 self.add_global_vectors("towards_image", "towards_object", "centre_normal", "support_normal", "majoraxis") 750 751 def _get_intersection(self, Ray): 752 """Return the intersection point between Ray and the ellipsoidal mirror surface.""" 753 ux, uy, uz = Ray.vector 754 xA, yA, zA = Ray.point 755 756 da = (uy**2 + uz**2) / self.b**2 + (ux / self.a) ** 2 757 db = 2 * ((uy * yA + uz * zA) / self.b**2 + (ux * xA) / self.a**2) 758 dc = (yA**2 + zA**2) / self.b**2 + (xA / self.a) ** 2 - 1 759 760 Solution = mgeo.SolverQuadratic(da, db, dc) 761 Solution = mgeo.KeepPositiveSolution(Solution) 762 763 ListPointIntersection = [] 764 C = self.get_centre_ref() 765 for t in Solution: 766 Intersect = Ray.vector * t + Ray.point 767 if Intersect[2] < 0 and Intersect - C in self.support: 768 ListPointIntersection.append(Intersect) 769 770 return _IntersectionRayMirror(Ray.point, ListPointIntersection) 771 772 def get_local_normal(self, Point): 773 """Return the normal unit vector on the ellipsoidal surface at point Point.""" 774 Gradient = Vector(np.zeros(3)) 775 776 Gradient[0] = -Point[0] / self.a**2 777 Gradient[1] = -Point[1] / self.b**2 778 Gradient[2] = -Point[2] / self.b**2 779 780 return Gradient.normalized() 781 782 def _get_centre_ref(self): 783 """Return 3D coordinates of the point on the mirror surface at the center of its support.""" 784 foci = 2 * np.sqrt(self.a**2 - self.b**2) # distance between the foci 785 h = -foci / 2 / np.tan(self._offaxisangle) 786 R = np.sqrt(foci**2 / 4 + h**2) 787 sign = 1 788 if math.isclose(self._offaxisangle, np.pi / 2): 789 h = 0 790 elif self._offaxisangle > np.pi / 2: 791 h = -h 792 sign = -1 793 a = 1 - self.a**2 / self.b**2 794 b = -2 * h 795 c = self.a**2 + h**2 - R**2 796 z = (-b + sign * np.sqrt(b**2 - 4 * a * c)) / (2 * a) 797 if math.isclose(z**2, self.b**2): 798 return np.array([0, 0, -self.b]) 799 x = self.a * np.sqrt(1 - z**2 / self.b**2) 800 centre = Point([x, 0, sign * z]) 801 return centre 802 803 def _zfunc(self, PointArray): 804 x = PointArray[:,0] 805 y = PointArray[:,1] 806 x-= self.r0[0] 807 y-= self.r0[1] 808 z = -np.sqrt(1 - (x / self.a)**2 - (y / self.b)**2) 809 z += self.r0[2] 810 return z 811 812 813# %% Cylindrical mirror definitions 814class MirrorCylindrical(Mirror): 815 """ 816 Cylindrical mirror surface with eqn. $y^2 + z^2 = R^2$, where $R$ is the radius. 817 818 Attributes 819 ---------- 820 radius : float 821 Radius of curvature. A postive value for concave mirror, a negative value for a convex mirror. 822 823 support : ART.ModuleSupport.Support 824 825 type : str 'Cylindrical Mirror'. 826 827 Methods 828 ------- 829 MirrorCylindrical.get_normal(Point) 830 831 MirrorCylindrical.get_centre() 832 833 MirrorCylindrical.get_grid3D(NbPoints) 834 """ 835 836 def __init__(self, Support, Radius, **kwargs): 837 """ 838 Construct a cylindrical mirror. 839 840 Parameters 841 ---------- 842 Radius : float 843 The radius of curvature in mm. A postive value for concave mirror, a negative value for a convex mirror. 844 845 Support : ART.ModuleSupport.Support 846 847 """ 848 super().__init__() 849 if Radius < 0: 850 self.type = "CylindricalCX Mirror" 851 self.curvature = Curvature.CONVEX 852 self.radius = -Radius 853 else: 854 self.type = "CylindricalCC Mirror" 855 self.curvature = Curvature.CONCAVE 856 self.radius = Radius 857 858 self.support = Support 859 860 if "Surface" in kwargs: 861 self.Surface = kwargs["Surface"] 862 863 self.r0 = Point([0.0, 0.0, -self.radius]) 864 865 self.support_normal_ref = Vector([0, 0, 1.0]) 866 self.majoraxis_ref = Vector([1.0, 0, 0]) 867 self.centre_ref = Point([0.0, 0.0, -self.radius]) 868 869 self.add_global_points("centre") 870 self.add_global_vectors("support_normal", "majoraxis") 871 872 def _get_intersection(self, Ray): 873 """Return the intersection point between the Ray and the cylinder.""" 874 uy = Ray.vector[1] 875 uz = Ray.vector[2] 876 yA = Ray.point[1] 877 zA = Ray.point[2] 878 879 a = uy**2 + uz**2 880 b = 2 * (uy * yA + uz * zA) 881 c = yA**2 + zA**2 - self.radius**2 882 883 Solution = mgeo.SolverQuadratic(a, b, c) 884 Solution = mgeo.KeepPositiveSolution(Solution) 885 886 ListPointIntersection = [] 887 for t in Solution: 888 Intersect = Ray.vector * t + Ray.point 889 if Intersect[2] < 0 and Intersect in self.support: 890 ListPointIntersection.append(Intersect) 891 892 return _IntersectionRayMirror(Ray.point, ListPointIntersection) 893 894 def get_local_normal(self, Point): 895 """Return the normal unit vector on the cylinder surface at point P.""" 896 Gradient = Vector([0, -Point[1], -Point[2]]) 897 return Gradient.normalized() 898 899 def get_centre(self): 900 """Return 3D coordinates of the point on the mirror surface at the center of its support.""" 901 return Point([0, 0, -self.radius]) 902 903 def _zfunc(self, PointArray): 904 y = PointArray[:,1] 905 return -np.sqrt(self.radius**2 - y**2) 906 907 908# %% Grazing parabolic mirror definitions 909# A grazing parabola is no different from a parabola, it's just a question of what we consider to be the support. 910# For an ordinary parabola, the support is the xy-plane, for a grazing parabola, the support is parallel to the surface at the optical center. 911 912class GrazingParabola(Mirror): 913 r""" 914 A parabolic mirror with a support parallel to the surface at the optical center. 915 Needs to be completed, currently not functional. 916 TODO 917 """ 918 919 def __init__(self, Support, 920 FocalEffective: float=None, 921 OffAxisAngle: float = None, 922 FocalParent: float = None, 923 RadiusParent: float = None, 924 OffAxisDistance: float = None, 925 MoreThan90: bool = None, 926 **kwargs): 927 """ 928 Initialise a Parabolic mirror. 929 930 Parameters 931 ---------- 932 FocalEffective : float 933 Effective focal length of the parabola in mm. 934 935 OffAxisAngle : float 936 Off-axis angle *in degrees* of the parabola. 937 938 Support : ART.ModuleSupport.Support 939 940 """ 941 super().__init__() 942 self.curvature = Curvature.CONCAVE 943 self.support = Support 944 self.type = "Grazing Parabolic Mirror" 945 946 if "Surface" in kwargs: 947 self.Surface = kwargs["Surface"] 948 949 parameter_calculator = OAP_calculator() 950 values,steps = parameter_calculator.calculate_values(verify_consistency=True, 951 fs=FocalEffective, 952 theta=OffAxisAngle, 953 fp=FocalParent, 954 Rc=RadiusParent, 955 OAD=OffAxisDistance, 956 more_than_90=MoreThan90) 957 self._offaxisangle = values["theta"] 958 self._feff = values["fs"] 959 self._p = values["p"] 960 self._offaxisdistance = values["OAD"] 961 self._fparent = values["fp"] 962 self._rparent = values["Rc"] 963 964 #self.r0 = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 965 self.r0 = Point([0, 0, 0]) 966 967 self.centre_ref = Point([0, 0, 0]) 968 self.support_normal_ref = Vector([0, 0, 1.0]) 969 self.majoraxis_ref = Vector([1.0, 0, 0]) 970 971 focus_ref = Point([np.sqrt(self._feff**2 - self._offaxisdistance**2), 0, self._offaxisdistance]) # frame where x is as collimated direction 972 # We need to rotate it in the trigonometric direction by the 90°-theta/2 973 angle = np.pi/2 - self._offaxisangle/2 974 self.focus_ref = Point([focus_ref[0]*np.cos(angle)+focus_ref[2]*np.sin(angle),focus_ref[1], -focus_ref[0]*np.sin(angle)+focus_ref[2]*np.cos(angle), ]) 975 976 self.towards_focusing_ref = (self.focus_ref-self.centre_ref).normalized() 977 towards_collimated_ref = Vector([-1.0, 0, 0]) # Again, we need to rotate it by the same angle 978 self.collimated_ref = Point([-np.cos(angle), 0, np.sin(angle)]) 979 980 self.add_global_points("focus", "centre") 981 self.add_global_vectors("towards_focusing", "towards_collimated", "support_normal", "majoraxis") 982 983 def _get_intersection(self, Ray): 984 """ 985 Return the intersection point between the ray and the parabola. 986 """ 987 pass 988 989 def get_local_normal(self, Point): 990 """ 991 Return the normal unit vector on the paraboloid surface at point Point. 992 """ 993 pass 994 995 def _zfunc(self, PointArray): 996 pass 997 998 999# %% Reflections on mirrors, is it useful to keep this? 1000def _ReflectionMirrorRay(Mirror, PointMirror, Ray): 1001 """ 1002 Return the reflected ray according to the law of reflection. 1003 1004 Parameters 1005 ---------- 1006 Mirror : Mirror-objectS 1007 1008 PointMirror : np.ndarray 1009 Point of reflection on the mirror surface. 1010 1011 Ray : Ray-object 1012 1013 """ 1014 PointRay = Ray.point 1015 VectorRay = Ray.vector 1016 NormalMirror = Mirror.get_local_normal(PointMirror) 1017 1018 #VectorRayReflected = mgeo.SymmetricalVector(-VectorRay, NormalMirror) 1019 VectorRayReflected = VectorRay- 2*NormalMirror*np.dot(VectorRay,NormalMirror) # Is it any better than SymmetricalVector? 1020 1021 RayReflected = Ray.copy_ray() 1022 RayReflected.point = PointMirror 1023 RayReflected.vector = VectorRayReflected 1024 RayReflected.incidence = mgeo.AngleBetweenTwoVectors( 1025 -VectorRay, NormalMirror 1026 ) 1027 RayReflected.path = Ray.path + (np.linalg.norm(PointMirror - PointRay),) 1028 1029 return RayReflected 1030 1031 1032def ReflectionMirrorRayList(Mirror, ListRay, IgnoreDefects=False): 1033 """ 1034 Return the the reflected rays according to the law of reflection for the list of incident rays ListRay. 1035 1036 Rays that do not hit the support are not further propagated. 1037 1038 Updates the reflected rays' incidence angle and path. 1039 1040 Parameters 1041 ---------- 1042 Mirror : Mirror-object 1043 1044 ListRay : list[Ray-object] 1045 1046 """ 1047 Deformed = type(Mirror) == DeformedMirror 1048 ListRayReflected = [] 1049 for k in ListRay: 1050 PointMirror = Mirror._get_intersection(k) 1051 1052 if PointMirror is not None: 1053 if Deformed and IgnoreDefects: 1054 M = Mirror.Mirror 1055 else: 1056 M = Mirror 1057 RayReflected = _ReflectionMirrorRay(M, PointMirror, k) 1058 ListRayReflected.append(RayReflected) 1059 return ListRayReflected 1060 1061 1062# %% Deformed mirror definitions, need to be reworked, maybe implemented as coatings? 1063class DeformedMirror(Mirror): 1064 def __init__(self, Mirror, DeformationList): 1065 self.Mirror = Mirror 1066 self.DeformationList = DeformationList 1067 self.type = Mirror.type 1068 self.support = self.Mirror.support 1069 1070 def get_local_normal(self, PointMirror): 1071 base_normal = self.Mirror.get_normal(PointMirror) 1072 C = self.get_centre() 1073 defects_normals = [ 1074 d.get_normal(PointMirror - C) for d in self.DeformationList 1075 ] 1076 for i in defects_normals: 1077 base_normal = mgeo.normal_add(base_normal, i) 1078 base_normal /= np.linalg.norm(base_normal) 1079 return base_normal 1080 1081 def get_centre(self): 1082 return self.Mirror.get_centre() 1083 1084 def get_grid3D(self, NbPoint, **kwargs): 1085 return self.Mirror.get_grid3D(NbPoint, **kwargs) 1086 1087 def _get_intersection(self, Ray): 1088 Intersect = self.Mirror._get_intersection(Ray) 1089 if Intersect is not None: 1090 h = sum( 1091 D.get_offset(Intersect - self.get_centre()) 1092 for D in self.DeformationList 1093 ) 1094 alpha = mgeo.AngleBetweenTwoVectors( 1095 -Ray.vector, self.Mirror.get_normal(Intersect) 1096 ) 1097 Intersect -= Ray.vector * h / np.cos(alpha) 1098 return Intersect
logger =
<Logger ARTcore.ModuleMirror (WARNING)>
class
Curvature(enum.Enum):
An enumeration.
CONVEX =
<Curvature.CONVEX: -1>
FLAT =
<Curvature.FLAT: 0>
CONCAVE =
<Curvature.CONCAVE: 1>
44class Mirror(moe.OpticalElement): 45 """Abstract base class for mirrors.""" 46 vectorised = False 47 def __init__(self): 48 self.Surface = msurf.IdealSurface() 49 super().__init__() 50 51 @abstractmethod 52 def _get_intersection(self, Ray): 53 """ 54 This method should return the intersection point between the ray and the mirror surface. 55 The calculation should be done in the reference frame of the mirror. 56 It essentially defines the mirror surface (but not its normal). 57 """ 58 pass 59 60 @abstractmethod 61 def get_local_normal(self, Point: np.ndarray): 62 """ 63 This method should return the normal unit vector in point 'Point' on the mirror surface. 64 The normal is in the reference frame of the mirror. 65 """ 66 pass 67 68 def propagate_raylist(self, RayList, alignment=False): 69 """ 70 Propagate a list of rays to the mirror and return the reflected rays. 71 72 Parameters 73 ---------- 74 RayList : a RayList (a class defined in ARTcore.ModuleOpticalRay) 75 Rays to propagate through the mirror. 76 77 alignment : bool 78 If True, the rays are not blocked by the support. Not implemented yet 79 80 Returns 81 ------- 82 ReflectedRayList : a RayList 83 The rays that were reflected by the mirror. 84 """ 85 #local_rays = [ray.to_basis(*self.basis) for ray in RayList] 86 local_rays = RayList.to_basis(*self.basis) 87 # So local_rays is a list of rays in the mirror reference frame 88 if self.vectorised and len(local_rays) > 10: 89 # If the class implements the vectorised version of the intersection and normal calculation 90 intersection_points, OK = self._get_intersections(local_rays) 91 # NB: the _get_intersections method should return two things: 92 # - the intersection points in the mirror reference frame with nan's for rays that do not intersect 93 # - a boolean array indicating which rays intersections make sense 94 local_normals = mgeo.VectorArray(np.zeros_like(intersection_points)) 95 local_normals[OK] = self.get_local_normals(intersection_points[OK]) 96 local_reflected_rays = [ 97 self.Surface.reflect_ray(local_rays[i], intersection_points[i], local_normals[i]) 98 for i in range(len(local_rays)) if OK[i]] 99 else: 100 intersection_points, OK = zip(*[self._get_intersection(r) for r in local_rays]) 101 # OK is a boolean list indicating if the ray intersects the mirror (for instance is in the support) 102 local_reflected_rays = [ 103 self.Surface.reflect_ray(local_rays[i], intersection_points[i], self.get_local_normal(intersection_points[i])) 104 for i in range(len(local_rays)) if OK[i]] 105 reflected_rays = mray.RayList.from_list(local_reflected_rays).from_basis(*self.basis) 106 if alignment and len(reflected_rays) == 0: 107 logger.error("Unexpected error occurred: no rays were reflected during alignment procedure.") 108 logger.error("This should not happen.") 109 logger.error(f"Mirror on which error ocurred: {self}") 110 logger.error(f"Alignment ray: {RayList[0]}") 111 logger.error(f"Alignment ray in mirror reference frame: {local_rays[0]}") 112 logger.error(f"Intersection point: {self._get_intersection(local_rays[0])[0]}") 113 logger.error(f"Support: {self.support}") 114 elif len(reflected_rays) == 0: 115 logger.warning("No rays were reflected by the mirror.") 116 logger.debug(f"Mirror: {self}") 117 logger.debug(f"First ray: {RayList[0]}") 118 logger.debug(f"First ray in mirror reference frame: {RayList[0].to_basis(self.r0, self.r, self.q)}") 119 logger.debug(f"First ray intersection point: {self._get_intersection(RayList[0].to_basis(self.r0, self.r, self.q))}") 120 return reflected_rays
Abstract base class for mirrors.
@abstractmethod
def
get_local_normal(self, Point: numpy.ndarray):
60 @abstractmethod 61 def get_local_normal(self, Point: np.ndarray): 62 """ 63 This method should return the normal unit vector in point 'Point' on the mirror surface. 64 The normal is in the reference frame of the mirror. 65 """ 66 pass
This method should return the normal unit vector in point 'Point' on the mirror surface. The normal is in the reference frame of the mirror.
def
propagate_raylist(self, RayList, alignment=False):
68 def propagate_raylist(self, RayList, alignment=False): 69 """ 70 Propagate a list of rays to the mirror and return the reflected rays. 71 72 Parameters 73 ---------- 74 RayList : a RayList (a class defined in ARTcore.ModuleOpticalRay) 75 Rays to propagate through the mirror. 76 77 alignment : bool 78 If True, the rays are not blocked by the support. Not implemented yet 79 80 Returns 81 ------- 82 ReflectedRayList : a RayList 83 The rays that were reflected by the mirror. 84 """ 85 #local_rays = [ray.to_basis(*self.basis) for ray in RayList] 86 local_rays = RayList.to_basis(*self.basis) 87 # So local_rays is a list of rays in the mirror reference frame 88 if self.vectorised and len(local_rays) > 10: 89 # If the class implements the vectorised version of the intersection and normal calculation 90 intersection_points, OK = self._get_intersections(local_rays) 91 # NB: the _get_intersections method should return two things: 92 # - the intersection points in the mirror reference frame with nan's for rays that do not intersect 93 # - a boolean array indicating which rays intersections make sense 94 local_normals = mgeo.VectorArray(np.zeros_like(intersection_points)) 95 local_normals[OK] = self.get_local_normals(intersection_points[OK]) 96 local_reflected_rays = [ 97 self.Surface.reflect_ray(local_rays[i], intersection_points[i], local_normals[i]) 98 for i in range(len(local_rays)) if OK[i]] 99 else: 100 intersection_points, OK = zip(*[self._get_intersection(r) for r in local_rays]) 101 # OK is a boolean list indicating if the ray intersects the mirror (for instance is in the support) 102 local_reflected_rays = [ 103 self.Surface.reflect_ray(local_rays[i], intersection_points[i], self.get_local_normal(intersection_points[i])) 104 for i in range(len(local_rays)) if OK[i]] 105 reflected_rays = mray.RayList.from_list(local_reflected_rays).from_basis(*self.basis) 106 if alignment and len(reflected_rays) == 0: 107 logger.error("Unexpected error occurred: no rays were reflected during alignment procedure.") 108 logger.error("This should not happen.") 109 logger.error(f"Mirror on which error ocurred: {self}") 110 logger.error(f"Alignment ray: {RayList[0]}") 111 logger.error(f"Alignment ray in mirror reference frame: {local_rays[0]}") 112 logger.error(f"Intersection point: {self._get_intersection(local_rays[0])[0]}") 113 logger.error(f"Support: {self.support}") 114 elif len(reflected_rays) == 0: 115 logger.warning("No rays were reflected by the mirror.") 116 logger.debug(f"Mirror: {self}") 117 logger.debug(f"First ray: {RayList[0]}") 118 logger.debug(f"First ray in mirror reference frame: {RayList[0].to_basis(self.r0, self.r, self.q)}") 119 logger.debug(f"First ray intersection point: {self._get_intersection(RayList[0].to_basis(self.r0, self.r, self.q))}") 120 return reflected_rays
Propagate a list of rays to the mirror and return the reflected rays.
Parameters
RayList : a RayList (a class defined in ARTcore.ModuleOpticalRay)
Rays to propagate through the mirror.
alignment : bool
If True, the rays are not blocked by the support. Not implemented yet
Returns
ReflectedRayList : a RayList
The rays that were reflected by the mirror.
def
getClosestSphere(Mirror, Npoints=1000):
445def GetClosestSphere(Mirror, Npoints=1000): 446 """ 447 This function calculates the closest sphere to the surface of a mirror. 448 It does so by sampling the surface of the mirror at Npoints points. 449 It then calculates the closest sphere to these points. 450 It returns the radius of the sphere and the center of the sphere. 451 """ 452 Points = mpm.sample_support(Mirror.support, Npoints=1000) 453 Points += Mirror.r0[:2] 454 Z = Mirror._zfunc(Points) 455 Points = mgeo.PointArray([Points[:, 0], Points[:, 1], Z]).T 456 spX, spY, spZ = Points[:, 0], Points[:, 1], Points[:, 2] 457 Center, Radius = BestFitSphere(spX, spY, spZ) 458 return Center, Radius
This function calculates the closest sphere to the surface of a mirror. It does so by sampling the surface of the mirror at Npoints points. It then calculates the closest sphere to these points. It returns the radius of the sphere and the center of the sphere.
def
getAsphericity(Mirror, Npoints=1000):
460def GetAsphericity(Mirror, Npoints=1000): 461 """ 462 This function calculates the maximum distance of the mirror surface to the closest sphere. 463 """ 464 center, radius = GetClosestSphere(Mirror, Npoints) 465 Points = mpm.sample_support(Mirror.support, Npoints=1000) 466 Points += Mirror.r0[:2] 467 Z = Mirror._zfunc(Points) 468 Points = mgeo.PointArray([Points[:, 0], Points[:, 1], Z]).T 469 Points_centered = Points - center 470 Distance = np.linalg.norm(Points_centered, axis=1) - radius 471 Distance*=1e3 # To convert to µm 472 return np.ptp(Distance)
This function calculates the maximum distance of the mirror surface to the closest sphere.
def
drawAsphericity(Mirror, Npoints=1000):
651def DrawAsphericity(Mirror, Npoints=1000): 652 """ 653 This function displays a map of the asphericity of the mirror. 654 It's a scatter plot of the points of the mirror surface, with the color representing the distance to the closest sphere. 655 The closest sphere is calculated by the function get_closest_sphere, so least square method. 656 657 Parameters 658 ---------- 659 Mirror : Mirror 660 The mirror to analyse. 661 662 Npoints : int, optional 663 The number of points to sample on the mirror surface. The default is 1000. 664 665 Returns 666 ------- 667 fig : Figure 668 The figure of the plot. 669 """ 670 plt.ion() 671 fig = plt.figure() 672 ax = Mirror.support._ContourSupport(fig) 673 center, radius = man.GetClosestSphere(Mirror, Npoints) 674 Points = mpm.sample_support(Mirror.support, Npoints=1000) 675 Points += Mirror.r0[:2] 676 Z = Mirror._zfunc(Points) 677 Points = mgeo.PointArray([Points[:, 0], Points[:, 1], Z]).T 678 X, Y = Points[:, 0] - Mirror.r0[0], Points[:, 1] - Mirror.r0[1] 679 Points_centered = Points - center 680 Distance = np.linalg.norm(Points_centered, axis=1) - radius 681 Distance*=1e3 # To convert to µm 682 p = plt.scatter(X, Y, c=Distance, s=15) 683 divider = man.make_axes_locatable(ax) 684 cax = divider.append_axes("right", size="5%", pad=0.05) 685 cbar = fig.colorbar(p, cax=cax) 686 cbar.set_label("Distance to closest sphere (µm)") 687 ax.set_xlabel("x (mm)") 688 ax.set_ylabel("y (mm)") 689 plt.title("Asphericity map", loc="right") 690 plt.tight_layout() 691 692 bbox = ax.get_position() 693 bbox.set_points(bbox.get_points() - np.array([[0.01, 0], [0.01, 0]])) 694 ax.set_position(bbox) 695 plt.show() 696 return fig
This function displays a map of the asphericity of the mirror. It's a scatter plot of the points of the mirror surface, with the color representing the distance to the closest sphere. The closest sphere is calculated by the function get_closest_sphere, so least square method.
Parameters
Mirror : Mirror The mirror to analyse.
Npoints : int, optional The number of points to sample on the mirror surface. The default is 1000.
Returns
fig : Figure The figure of the plot.
138class MirrorPlane(Mirror): 139 """ 140 A plane mirror. 141 142 Attributes 143 ---------- 144 support : ART.ModuleSupport.Support 145 Mirror support 146 type : str = 'Plane Mirror' 147 Human readable mirror type 148 Methods 149 ------- 150 MirrorPlane.get_normal(Point) 151 152 MirrorPlane.get_centre() 153 154 MirrorPlane.get_grid3D(NbPoints) 155 """ 156 vectorised = True 157 def __init__(self, Support, **kwargs): 158 """ 159 Create a plane mirror on a given Support. 160 161 Parameters 162 ---------- 163 Support : ART.ModuleSupport.Support 164 165 """ 166 super().__init__() 167 self.support = Support 168 self.type = "Plane Mirror" 169 self.curvature = Curvature.FLAT 170 171 if "Surface" in kwargs: 172 self.Surface = kwargs["Surface"] 173 174 self.r0 = Point([0.0, 0.0, 0.0]) 175 176 self.centre_ref = Point([0.0, 0.0, 0.0]) 177 self.support_normal_ref = Vector([0, 0, 1.0]) 178 self.majoraxis_ref = Vector([1.0, 0, 0]) 179 180 self.add_global_points("centre") 181 self.add_global_vectors("support_normal", "majoraxis") 182 183 def _get_intersection(self, Ray): 184 """Return the intersection point between Ray and the xy-plane.""" 185 t = -Ray.point[2] / Ray.vector[2] 186 PointIntersection = Ray.point + Ray.vector * t 187 return PointIntersection, t > 0 and PointIntersection in self.support 188 189 def get_local_normal(self, Point: np.ndarray): 190 """Return normal unit vector in point 'Point' on the plane mirror.""" 191 return Vector([0, 0, 1]) 192 193 def _get_intersections(self, RayList): 194 """ 195 Vectorised version of the intersection calculation. 196 """ 197 vectors = np.array([ray.vector for ray in RayList]) 198 points = np.array([ray.point for ray in RayList]) 199 200 t = -points[:,2] / vectors[:,2] 201 Points = points + vectors * t[:,np.newaxis] 202 OK = (t > 0) & np.array([p in self.support for p in Points]) 203 return mgeo.PointArray(Points), OK 204 205 def get_local_normals(self, Points): 206 return mgeo.VectorArray(np.zeros_like(Points) + [0, 0, 1]) 207 208 def _zfunc(self, PointArray): 209 x = PointArray[:,0] 210 y = PointArray[:,1] 211 return np.zeros_like(x)
A plane mirror.
Attributes
support : ART.ModuleSupport.Support
Mirror support
type : str = 'Plane Mirror'
Human readable mirror type
Methods
MirrorPlane.get_normal(Point)
MirrorPlane.get_centre()
MirrorPlane.get_grid3D(NbPoints)
MirrorPlane(Support, **kwargs)
157 def __init__(self, Support, **kwargs): 158 """ 159 Create a plane mirror on a given Support. 160 161 Parameters 162 ---------- 163 Support : ART.ModuleSupport.Support 164 165 """ 166 super().__init__() 167 self.support = Support 168 self.type = "Plane Mirror" 169 self.curvature = Curvature.FLAT 170 171 if "Surface" in kwargs: 172 self.Surface = kwargs["Surface"] 173 174 self.r0 = Point([0.0, 0.0, 0.0]) 175 176 self.centre_ref = Point([0.0, 0.0, 0.0]) 177 self.support_normal_ref = Vector([0, 0, 1.0]) 178 self.majoraxis_ref = Vector([1.0, 0, 0]) 179 180 self.add_global_points("centre") 181 self.add_global_vectors("support_normal", "majoraxis")
Create a plane mirror on a given Support.
Parameters
Support : ART.ModuleSupport.Support
def
get_local_normal(self, Point: numpy.ndarray):
189 def get_local_normal(self, Point: np.ndarray): 190 """Return normal unit vector in point 'Point' on the plane mirror.""" 191 return Vector([0, 0, 1])
Return normal unit vector in point 'Point' on the plane mirror.
def
render( Mirror, Npoints=1000, draw_support=False, draw_points=False, draw_vectors=True, recenter_support=True):
234def _RenderMirror(Mirror, Npoints=1000, draw_support = False, draw_points = False, draw_vectors = True , recenter_support = True): 235 """ 236 This function renders a mirror in 3D. 237 """ 238 mesh = Mirror._Render(Npoints) 239 p = pv.Plotter() 240 p.add_mesh(mesh) 241 if draw_support: 242 support = Mirror.support._Render() 243 if recenter_support: 244 support.translate(Mirror.r0,inplace=True) 245 p.add_mesh(support, color="gray", opacity=0.5) 246 247 if draw_vectors: 248 # We draw the important vectors of the optical element 249 # For that, if we have a "vectors" attribute, we use that 250 # (a dictionary with the vector names as keys and the colors as values) 251 # Otherwise we use the default vectors: "support_normal", "majoraxis" 252 if hasattr(Mirror, "vectors"): 253 vectors = Mirror.vectors 254 else: 255 vectors = {"support_normal_ref": "red", "majoraxis_ref": "blue"} 256 for vector, color in vectors.items(): 257 if hasattr(Mirror, vector): 258 p.add_arrows(Mirror.r0, 10*getattr(Mirror, vector), color=color) 259 260 if draw_points: 261 # We draw the important points of the optical element 262 # For that, if we have a "points" attribute, we use that 263 # (a dictionary with the point names as keys and the colors as values) 264 # Otherwise we use the default points: "centre_ref" 265 if hasattr(Mirror, "points"): 266 points = Mirror.points 267 else: 268 points = {"centre_ref": "red"} 269 for point, color in points.items(): 270 if hasattr(Mirror, point): 271 p.add_mesh(pv.Sphere(radius=1, center=getattr(Mirror, point)), color=color) 272 else: 273 p.add_mesh(pv.Sphere(radius=1, center=Mirror.r0), color="red") 274 p.show() 275 return p
This function renders a mirror in 3D.
Inherited Members
214class MirrorSpherical(Mirror): 215 """ 216 Spherical mirror surface with eqn. $x^2 + y^2 + z^2 = R^2$, where $R$ is the radius. 217 218  219 220 Attributes 221 ---------- 222 radius : float 223 Radius of curvature. A postive value for concave mirror, a negative value for a convex mirror. 224 225 support : ART.ModuleSupport.Support 226 227 type : str SphericalCC Mirror' or SphericalCX Mirror' 228 229 Methods 230 ------- 231 MirrorSpherical.get_normal(Point) 232 233 MirrorSpherical.get_centre() 234 235 MirrorSpherical.get_grid3D(NbPoints) 236 237 """ 238 def __init__(self, Support, Radius, **kwargs): 239 """ 240 Construct a spherical mirror. 241 242 Parameters 243 ---------- 244 Radius : float 245 The radius of curvature in mm. A postive value for concave mirror, a negative value for a convex mirror. 246 247 Support : ART.ModuleSupport.Support 248 249 """ 250 super().__init__() 251 if Radius < 0: 252 self.type = "SphericalCX Mirror" 253 self.curvature = Curvature.CONVEX 254 self.radius = -Radius 255 elif Radius > 0: 256 self.type = "SphericalCC Mirror" 257 self.curvature = Curvature.CONCAVE 258 self.radius = Radius 259 else: 260 raise ValueError("Radius of curvature must be non-zero") 261 self.support = Support 262 263 if "Surface" in kwargs: 264 self.Surface = kwargs["Surface"] 265 266 self.r0 = Point([0.0, 0, -self.radius]) 267 268 self.centre_ref = Point([0.0, 0, -self.radius]) 269 self.sphere_center_ref = Point([0.0, 0, 0]) 270 self.focus_ref = Point([0.0, 0, -self.radius/2]) 271 272 self.support_normal_ref = Vector([0, 0, 1.0]) 273 self.majoraxis_ref = Vector([1.0, 0, 0]) 274 self.towards_focus_ref = (self.focus_ref - self.centre_ref).normalized() 275 276 self.add_global_points("sphere_center", "focus", "centre") 277 self.add_global_vectors("towards_focus", "majoraxis", "support_normal") 278 279 def _get_intersection(self, Ray): 280 """Return the intersection point between the ray and the sphere.""" 281 a = np.dot(Ray.vector, Ray.vector) 282 b = 2 * np.dot(Ray.vector, Ray.point) 283 c = np.dot(Ray.point, Ray.point) - self.radius**2 284 285 Solution = mgeo.SolverQuadratic(a, b, c) 286 Solution = mgeo.KeepPositiveSolution(Solution) 287 288 ListPointIntersection = [] 289 for t in Solution: 290 Intersect = Ray.vector * t + Ray.point 291 if Intersect[2] < 0 and Intersect in self.support: 292 ListPointIntersection.append(Intersect) 293 294 return _IntersectionRayMirror(Ray.point, ListPointIntersection) 295 296 def get_local_normal(self, Point): 297 """Return the normal unit vector on the spherical surface at point Point.""" 298 return -Vector(Point).normalized() 299 300 def _zfunc(self, PointArray): 301 x = PointArray[:,0] 302 y = PointArray[:,1] 303 return -np.sqrt(self.radius**2 - x**2 - y**2)
Spherical mirror surface with eqn. $x^2 + y^2 + z^2 = R^2$, where $R$ is the radius.
Attributes
radius : float
Radius of curvature. A postive value for concave mirror, a negative value for a convex mirror.
support : ART.ModuleSupport.Support
type : str SphericalCC Mirror' or SphericalCX Mirror'
Methods
MirrorSpherical.get_normal(Point)
MirrorSpherical.get_centre()
MirrorSpherical.get_grid3D(NbPoints)
MirrorSpherical(Support, Radius, **kwargs)
238 def __init__(self, Support, Radius, **kwargs): 239 """ 240 Construct a spherical mirror. 241 242 Parameters 243 ---------- 244 Radius : float 245 The radius of curvature in mm. A postive value for concave mirror, a negative value for a convex mirror. 246 247 Support : ART.ModuleSupport.Support 248 249 """ 250 super().__init__() 251 if Radius < 0: 252 self.type = "SphericalCX Mirror" 253 self.curvature = Curvature.CONVEX 254 self.radius = -Radius 255 elif Radius > 0: 256 self.type = "SphericalCC Mirror" 257 self.curvature = Curvature.CONCAVE 258 self.radius = Radius 259 else: 260 raise ValueError("Radius of curvature must be non-zero") 261 self.support = Support 262 263 if "Surface" in kwargs: 264 self.Surface = kwargs["Surface"] 265 266 self.r0 = Point([0.0, 0, -self.radius]) 267 268 self.centre_ref = Point([0.0, 0, -self.radius]) 269 self.sphere_center_ref = Point([0.0, 0, 0]) 270 self.focus_ref = Point([0.0, 0, -self.radius/2]) 271 272 self.support_normal_ref = Vector([0, 0, 1.0]) 273 self.majoraxis_ref = Vector([1.0, 0, 0]) 274 self.towards_focus_ref = (self.focus_ref - self.centre_ref).normalized() 275 276 self.add_global_points("sphere_center", "focus", "centre") 277 self.add_global_vectors("towards_focus", "majoraxis", "support_normal")
Construct a spherical mirror.
Parameters
Radius : float
The radius of curvature in mm. A postive value for concave mirror, a negative value for a convex mirror.
Support : ART.ModuleSupport.Support
def
get_local_normal(self, Point):
296 def get_local_normal(self, Point): 297 """Return the normal unit vector on the spherical surface at point Point.""" 298 return -Vector(Point).normalized()
Return the normal unit vector on the spherical surface at point Point.
def
render( Mirror, Npoints=1000, draw_support=False, draw_points=False, draw_vectors=True, recenter_support=True):
234def _RenderMirror(Mirror, Npoints=1000, draw_support = False, draw_points = False, draw_vectors = True , recenter_support = True): 235 """ 236 This function renders a mirror in 3D. 237 """ 238 mesh = Mirror._Render(Npoints) 239 p = pv.Plotter() 240 p.add_mesh(mesh) 241 if draw_support: 242 support = Mirror.support._Render() 243 if recenter_support: 244 support.translate(Mirror.r0,inplace=True) 245 p.add_mesh(support, color="gray", opacity=0.5) 246 247 if draw_vectors: 248 # We draw the important vectors of the optical element 249 # For that, if we have a "vectors" attribute, we use that 250 # (a dictionary with the vector names as keys and the colors as values) 251 # Otherwise we use the default vectors: "support_normal", "majoraxis" 252 if hasattr(Mirror, "vectors"): 253 vectors = Mirror.vectors 254 else: 255 vectors = {"support_normal_ref": "red", "majoraxis_ref": "blue"} 256 for vector, color in vectors.items(): 257 if hasattr(Mirror, vector): 258 p.add_arrows(Mirror.r0, 10*getattr(Mirror, vector), color=color) 259 260 if draw_points: 261 # We draw the important points of the optical element 262 # For that, if we have a "points" attribute, we use that 263 # (a dictionary with the point names as keys and the colors as values) 264 # Otherwise we use the default points: "centre_ref" 265 if hasattr(Mirror, "points"): 266 points = Mirror.points 267 else: 268 points = {"centre_ref": "red"} 269 for point, color in points.items(): 270 if hasattr(Mirror, point): 271 p.add_mesh(pv.Sphere(radius=1, center=getattr(Mirror, point)), color=color) 272 else: 273 p.add_mesh(pv.Sphere(radius=1, center=Mirror.r0), color="red") 274 p.show() 275 return p
This function renders a mirror in 3D.
Inherited Members
306class MirrorParabolic(Mirror): 307 r""" 308 A paraboloid with vertex at the origin $O=[0,0,0]$ and symmetry axis z: 309 $z = \frac{1}{4f}[x^2 + y^2]$ where $f$ is the focal lenght of the *mother* 310 parabola (i.e. measured from its center at $O$ to the focal point $F$). 311 312 The center of the support is shifted along the x-direction by the off-axis distance $x_c$. 313 This leads to an *effective focal length* $f_\mathrm{eff}$, measured from the shifted center 314 of the support $P$ to the focal point $F$. 315 It is related to the mother focal length by $f = f_\\mathrm{eff} \cos^2(\alpha/2) $, 316 or equivalently $ p = 2f = f_\mathrm{eff} (1+\\cos\\alpha)$, where $\\alpha$ 317 is the off-axis angle, and $p = 2f$ is called the semi latus rectum. 318 319 Another useful relationship is that between the off-axis distance and the resulting 320 off-axis angle: $x_c = 2 f \tan(\alpha/2)$. 321 322 323  324 325 Attributes 326 ---------- 327 offaxisangle : float 328 Off-axis angle of the parabola. Modifying this also updates p, keeping feff constant. 329 Attention: The off-axis angle must be *given in degrees*, but is stored and will be *returned in radian* ! 330 331 feff : float 332 Effective focal length of the parabola in mm. Modifying this also updates p, keeping offaxisangle constant. 333 334 p : float 335 Semi latus rectum of the parabola in mm. Modifying this also updates feff, keeping offaxisangle constant. 336 337 support : ART.ModuleSupport.Support 338 339 type : str 'Parabolic Mirror'. 340 341 Methods 342 ------- 343 MirrorParabolic.get_normal(Point) 344 345 MirrorParabolic.get_centre() 346 347 MirrorParabolic.get_grid3D(NbPoints) 348 349 """ 350 #vectorised = True # Unitl I fix the closest solution thing 351 def __init__(self, Support, 352 FocalEffective: float=None, 353 OffAxisAngle: float = None, 354 FocalParent: float = None, 355 RadiusParent: float = None, 356 OffAxisDistance: float = None, 357 MoreThan90: bool = None, 358 **kwargs): 359 """ 360 Initialise a Parabolic mirror. 361 362 Parameters 363 ---------- 364 FocalEffective : float 365 Effective focal length of the parabola in mm. 366 367 OffAxisAngle : float 368 Off-axis angle *in degrees* of the parabola. 369 370 Support : ART.ModuleSupport.Support 371 372 """ 373 super().__init__() 374 self.curvature = Curvature.CONCAVE 375 self.support = Support 376 self.type = "Parabolic Mirror" 377 378 if "Surface" in kwargs: 379 self.Surface = kwargs["Surface"] 380 381 parameter_calculator = OAP_calculator() 382 values,steps = parameter_calculator.calculate_values(verify_consistency=True, 383 fs=FocalEffective, 384 theta=OffAxisAngle, 385 fp=FocalParent, 386 Rc=RadiusParent, 387 OAD=OffAxisDistance, 388 more_than_90=MoreThan90) 389 self._offaxisangle = values["theta"] 390 self._feff = values["fs"] 391 self._p = values["p"] 392 self._offaxisdistance = values["OAD"] 393 self._fparent = values["fp"] 394 self._rparent = values["Rc"] 395 396 self.r0 = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 397 398 self.centre_ref = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 399 self.support_normal_ref = Vector([0, 0, 1.0]) 400 self.majoraxis_ref = Vector([1.0, 0, 0]) 401 402 self.focus_ref = Point([0.0, 0, self._fparent]) 403 self.towards_focusing_ref = (self.focus_ref-self.centre_ref).normalized() 404 self.towards_collimated_ref = Vector([0, 0, 1.0]) 405 406 self.add_global_points("focus", "centre") 407 self.add_global_vectors("towards_focusing", "towards_collimated", "support_normal", "majoraxis") 408 409 def _get_intersection(self, Ray): 410 """Return the intersection point between the ray and the parabola.""" 411 ux, uy, uz = Ray.vector 412 xA, yA, zA = Ray.point 413 414 da = ux**2 + uy**2 415 db = 2 * (ux * xA + uy * yA) - 2 * self._p * uz 416 dc = xA**2 + yA**2 - 2 * self._p * zA 417 418 Solution = mgeo.SolverQuadratic(da, db, dc) 419 Solution = mgeo.KeepPositiveSolution(Solution) 420 IntersectionPoint = [Ray.point + Ray.vector * t for t in Solution if t > 0] 421 distances = [mgeo.Vector(i - self.r0).norm for i in IntersectionPoint] 422 IntersectionPoint = IntersectionPoint[distances.index(min(distances))] 423 return IntersectionPoint, IntersectionPoint-self.r0 in self.support 424 425 def _get_intersections(self, RayList): 426 """ 427 Vectorised version of the intersection calculation. 428 """ 429 vectors = np.array([ray.vector for ray in RayList]) 430 points = np.array([ray.point for ray in RayList]) 431 Solutions = np.zeros(len(RayList), dtype=float) 432 433 da = np.sum(vectors[:,0:2]**2, axis=1) 434 db = 2 * (np.sum(vectors[:,0:2] * points[:,0:2], axis=1) - self._p * vectors[:,2]) 435 dc = np.sum(points[:,0:2]**2, axis=1) - 2 * self._p * points[:,2] 436 437 linear = da == 0 438 Solutions[linear] = -dc[linear] / db[linear] 439 440 centered = (dc==0) & (da!=0) # If equation is of form ax**2 + bx = 0 441 Solutions[centered] = np.maximum(-db[centered] / da[centered], 0) 442 443 unstable = (4*da*dc) < db**2 * 1e-10 # Cases leading to catastrophic cancellation 444 445 Deltas = db**2 - 4 * da * dc 446 OK = Deltas >= 0 447 SolutionsPlus = (-db + np.sqrt(Deltas)) / 2 / da 448 SolutionsMinus = (-db - np.sqrt(Deltas)) / 2 / da 449 SolutionsPlus[unstable] = dc[unstable] / (da[unstable] * SolutionsMinus[unstable]) 450 SolutionsPlus = np.where(SolutionsPlus >= 0, SolutionsPlus, np.inf) 451 SolutionsMinus = np.where(SolutionsMinus >= 0, SolutionsMinus, np.inf) 452 Solutions = np.minimum(SolutionsPlus, SolutionsMinus) 453 Solutions = np.maximum(Solutions, 0) 454 Points = points + vectors * Solutions[:,np.newaxis] 455 OK = OK & (Solutions > 0) & np.array([p-self.r0 in self.support for p in Points]) 456 return mgeo.PointArray(Points), OK 457 458 def get_local_normal(self, Point): 459 """Return the normal unit vector on the paraboloid surface at point Point.""" 460 Gradient = Vector(np.zeros(3)) 461 Gradient[0] = -Point[0] 462 Gradient[1] = -Point[1] 463 Gradient[2] = self._p 464 return Gradient.normalized() 465 466 def get_local_normals(self, Points): 467 """ 468 Vectorised version of the normal calculation. 469 """ 470 Gradients = mgeo.VectorArray(np.zeros_like(Points)) 471 Gradients[:,0] = -Points[:,0] 472 Gradients[:,1] = -Points[:,1] 473 Gradients[:,2] = self._p 474 return Gradients.normalized() 475 476 def _zfunc(self, PointArray): 477 x = PointArray[:,0] 478 y = PointArray[:,1] 479 return (x**2 + y**2) / 2 / self._p
A paraboloid with vertex at the origin $O=[0,0,0]$ and symmetry axis z: $z = \frac{1}{4f}[x^2 + y^2]$ where $f$ is the focal lenght of the mother parabola (i.e. measured from its center at $O$ to the focal point $F$).
The center of the support is shifted along the x-direction by the off-axis distance $x_c$. This leads to an *effective focal length* $f_\mathrm{eff}$, measured from the shifted center of the support $P$ to the focal point $F$. It is related to the mother focal length by $f = f_\mathrm{eff} \cos^2(\alpha/2) $, or equivalently $ p = 2f = f_\mathrm{eff} (1+\cos\alpha)$, where $\alpha$ is the off-axis angle, and $p = 2f$ is called the semi latus rectum.
Another useful relationship is that between the off-axis distance and the resulting off-axis angle: $x_c = 2 f \tan(\alpha/2)$.
Attributes
offaxisangle : float
Off-axis angle of the parabola. Modifying this also updates p, keeping feff constant.
Attention: The off-axis angle must be *given in degrees*, but is stored and will be *returned in radian* !
feff : float
Effective focal length of the parabola in mm. Modifying this also updates p, keeping offaxisangle constant.
p : float
Semi latus rectum of the parabola in mm. Modifying this also updates feff, keeping offaxisangle constant.
support : ART.ModuleSupport.Support
type : str 'Parabolic Mirror'.
Methods
MirrorParabolic.get_normal(Point)
MirrorParabolic.get_centre()
MirrorParabolic.get_grid3D(NbPoints)
MirrorParabolic( Support, FocalEffective: float = None, OffAxisAngle: float = None, FocalParent: float = None, RadiusParent: float = None, OffAxisDistance: float = None, MoreThan90: bool = None, **kwargs)
351 def __init__(self, Support, 352 FocalEffective: float=None, 353 OffAxisAngle: float = None, 354 FocalParent: float = None, 355 RadiusParent: float = None, 356 OffAxisDistance: float = None, 357 MoreThan90: bool = None, 358 **kwargs): 359 """ 360 Initialise a Parabolic mirror. 361 362 Parameters 363 ---------- 364 FocalEffective : float 365 Effective focal length of the parabola in mm. 366 367 OffAxisAngle : float 368 Off-axis angle *in degrees* of the parabola. 369 370 Support : ART.ModuleSupport.Support 371 372 """ 373 super().__init__() 374 self.curvature = Curvature.CONCAVE 375 self.support = Support 376 self.type = "Parabolic Mirror" 377 378 if "Surface" in kwargs: 379 self.Surface = kwargs["Surface"] 380 381 parameter_calculator = OAP_calculator() 382 values,steps = parameter_calculator.calculate_values(verify_consistency=True, 383 fs=FocalEffective, 384 theta=OffAxisAngle, 385 fp=FocalParent, 386 Rc=RadiusParent, 387 OAD=OffAxisDistance, 388 more_than_90=MoreThan90) 389 self._offaxisangle = values["theta"] 390 self._feff = values["fs"] 391 self._p = values["p"] 392 self._offaxisdistance = values["OAD"] 393 self._fparent = values["fp"] 394 self._rparent = values["Rc"] 395 396 self.r0 = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 397 398 self.centre_ref = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 399 self.support_normal_ref = Vector([0, 0, 1.0]) 400 self.majoraxis_ref = Vector([1.0, 0, 0]) 401 402 self.focus_ref = Point([0.0, 0, self._fparent]) 403 self.towards_focusing_ref = (self.focus_ref-self.centre_ref).normalized() 404 self.towards_collimated_ref = Vector([0, 0, 1.0]) 405 406 self.add_global_points("focus", "centre") 407 self.add_global_vectors("towards_focusing", "towards_collimated", "support_normal", "majoraxis")
Initialise a Parabolic mirror.
Parameters
FocalEffective : float
Effective focal length of the parabola in mm.
OffAxisAngle : float
Off-axis angle *in degrees* of the parabola.
Support : ART.ModuleSupport.Support
def
get_local_normal(self, Point):
458 def get_local_normal(self, Point): 459 """Return the normal unit vector on the paraboloid surface at point Point.""" 460 Gradient = Vector(np.zeros(3)) 461 Gradient[0] = -Point[0] 462 Gradient[1] = -Point[1] 463 Gradient[2] = self._p 464 return Gradient.normalized()
Return the normal unit vector on the paraboloid surface at point Point.
def
get_local_normals(self, Points):
466 def get_local_normals(self, Points): 467 """ 468 Vectorised version of the normal calculation. 469 """ 470 Gradients = mgeo.VectorArray(np.zeros_like(Points)) 471 Gradients[:,0] = -Points[:,0] 472 Gradients[:,1] = -Points[:,1] 473 Gradients[:,2] = self._p 474 return Gradients.normalized()
Vectorised version of the normal calculation.
def
render( Mirror, Npoints=1000, draw_support=False, draw_points=False, draw_vectors=True, recenter_support=True):
234def _RenderMirror(Mirror, Npoints=1000, draw_support = False, draw_points = False, draw_vectors = True , recenter_support = True): 235 """ 236 This function renders a mirror in 3D. 237 """ 238 mesh = Mirror._Render(Npoints) 239 p = pv.Plotter() 240 p.add_mesh(mesh) 241 if draw_support: 242 support = Mirror.support._Render() 243 if recenter_support: 244 support.translate(Mirror.r0,inplace=True) 245 p.add_mesh(support, color="gray", opacity=0.5) 246 247 if draw_vectors: 248 # We draw the important vectors of the optical element 249 # For that, if we have a "vectors" attribute, we use that 250 # (a dictionary with the vector names as keys and the colors as values) 251 # Otherwise we use the default vectors: "support_normal", "majoraxis" 252 if hasattr(Mirror, "vectors"): 253 vectors = Mirror.vectors 254 else: 255 vectors = {"support_normal_ref": "red", "majoraxis_ref": "blue"} 256 for vector, color in vectors.items(): 257 if hasattr(Mirror, vector): 258 p.add_arrows(Mirror.r0, 10*getattr(Mirror, vector), color=color) 259 260 if draw_points: 261 # We draw the important points of the optical element 262 # For that, if we have a "points" attribute, we use that 263 # (a dictionary with the point names as keys and the colors as values) 264 # Otherwise we use the default points: "centre_ref" 265 if hasattr(Mirror, "points"): 266 points = Mirror.points 267 else: 268 points = {"centre_ref": "red"} 269 for point, color in points.items(): 270 if hasattr(Mirror, point): 271 p.add_mesh(pv.Sphere(radius=1, center=getattr(Mirror, point)), color=color) 272 else: 273 p.add_mesh(pv.Sphere(radius=1, center=Mirror.r0), color="red") 274 p.show() 275 return p
This function renders a mirror in 3D.
Inherited Members
483class MirrorToroidal(Mirror): 484 r""" 485 Toroidal mirror surface with eqn.$(\sqrt{x^2+z^2}-R)^2 + y^2 = r^2$ where $R$ and $r$ the major and minor radii. 486 487  488 489 Attributes 490 ---------- 491 majorradius : float 492 Major radius of the toroid in mm. 493 494 minorradius : float 495 Minor radius of the toroid in mm. 496 497 support : ART.ModuleSupport.Support 498 499 type : str 'Toroidal Mirror'. 500 501 Methods 502 ------- 503 MirrorToroidal.get_normal(Point) 504 505 MirrorToroidal.get_centre() 506 507 MirrorToroidal.get_grid3D(NbPoints) 508 509 """ 510 511 def __init__(self, Support, MajorRadius, MinorRadius, **kwargs): 512 """ 513 Construct a toroidal mirror. 514 515 Parameters 516 ---------- 517 MajorRadius : float 518 Major radius of the toroid in mm. 519 520 MinorRadius : float 521 Minor radius of the toroid in mm. 522 523 Support : ART.ModuleSupport.Support 524 525 """ 526 super().__init__() 527 self.support = Support 528 self.type = "Toroidal Mirror" 529 530 if "Surface" in kwargs: 531 self.Surface = kwargs["Surface"] 532 533 self.curvature = Curvature.CONCAVE 534 535 self.majorradius = MajorRadius 536 self.minorradius = MinorRadius 537 538 self.r0 = Point([0.0, 0.0, -MajorRadius - MinorRadius]) 539 540 self.centre_ref = Point([0.0, 0.0, -MajorRadius - MinorRadius]) 541 self.support_normal_ref = Vector([0, 0, 1.0]) 542 self.majoraxis_ref = Vector([1.0, 0, 0]) 543 544 self.add_global_vectors("support_normal", "majoraxis") 545 self.add_global_points("centre") 546 547 def _get_intersection(self, Ray): 548 """Return the intersection point between Ray and the toroidal mirror surface. The intersection points are given in the reference frame of the mirror""" 549 ux = Ray.vector[0] 550 uz = Ray.vector[2] 551 xA = Ray.point[0] 552 zA = Ray.point[2] 553 554 G = 4.0 * self.majorradius**2 * (ux**2 + uz**2) 555 H = 8.0 * self.majorradius**2 * (ux * xA + uz * zA) 556 I = 4.0 * self.majorradius**2 * (xA**2 + zA**2) 557 J = np.dot(Ray.vector, Ray.vector) 558 K = 2.0 * np.dot(Ray.vector, Ray.point) 559 L = ( 560 np.dot(Ray.point, Ray.point) 561 + self.majorradius**2 562 - self.minorradius**2 563 ) 564 565 a = J**2 566 b = 2 * J * K 567 c = 2 * J * L + K**2 - G 568 d = 2 * K * L - H 569 e = L**2 - I 570 571 Solution = mgeo.SolverQuartic(a, b, c, d, e) 572 if len(Solution) == 0: 573 return None, False 574 Solution = [t for t in Solution if t > 0] 575 IntersectionPoint = [Ray.point + Ray.vector *i for i in Solution] 576 distances = [mgeo.Vector(i - self.r0).norm for i in IntersectionPoint] 577 IntersectionPoint = IntersectionPoint[distances.index(min(distances))] 578 OK = IntersectionPoint - self.r0 in self.support and np.abs(IntersectionPoint[2]-self.r0[2]) < self.majorradius 579 return IntersectionPoint, OK 580 581 def get_local_normal(self, Point): 582 """Return the normal unit vector on the toroidal surface at point Point in the reference frame of the mirror""" 583 x = Point[0] 584 y = Point[1] 585 z = Point[2] 586 A = self.majorradius**2 - self.minorradius**2 587 588 Gradient = Vector(np.zeros(3)) 589 Gradient[0] = ( 590 4 * (x**3 + x * y**2 + x * z**2 + x * A) 591 - 8 * x * self.majorradius**2 592 ) 593 Gradient[1] = 4 * (y**3 + y * x**2 + y * z**2 + y * A) 594 Gradient[2] = ( 595 4 * (z**3 + z * x**2 + z * y**2 + z * A) 596 - 8 * z * self.majorradius**2 597 ) 598 599 return -Gradient.normalized() 600 601 def _zfunc(self, PointArray): 602 x = PointArray[:,0] 603 y = PointArray[:,1] 604 return -np.sqrt( 605 (np.sqrt(self.minorradius**2 - y**2) + self.majorradius)**2 - x**2 606 )
Toroidal mirror surface with eqn.$(\sqrt{x^2+z^2}-R)^2 + y^2 = r^2$ where $R$ and $r$ the major and minor radii.
Attributes
majorradius : float
Major radius of the toroid in mm.
minorradius : float
Minor radius of the toroid in mm.
support : ART.ModuleSupport.Support
type : str 'Toroidal Mirror'.
Methods
MirrorToroidal.get_normal(Point)
MirrorToroidal.get_centre()
MirrorToroidal.get_grid3D(NbPoints)
MirrorToroidal(Support, MajorRadius, MinorRadius, **kwargs)
511 def __init__(self, Support, MajorRadius, MinorRadius, **kwargs): 512 """ 513 Construct a toroidal mirror. 514 515 Parameters 516 ---------- 517 MajorRadius : float 518 Major radius of the toroid in mm. 519 520 MinorRadius : float 521 Minor radius of the toroid in mm. 522 523 Support : ART.ModuleSupport.Support 524 525 """ 526 super().__init__() 527 self.support = Support 528 self.type = "Toroidal Mirror" 529 530 if "Surface" in kwargs: 531 self.Surface = kwargs["Surface"] 532 533 self.curvature = Curvature.CONCAVE 534 535 self.majorradius = MajorRadius 536 self.minorradius = MinorRadius 537 538 self.r0 = Point([0.0, 0.0, -MajorRadius - MinorRadius]) 539 540 self.centre_ref = Point([0.0, 0.0, -MajorRadius - MinorRadius]) 541 self.support_normal_ref = Vector([0, 0, 1.0]) 542 self.majoraxis_ref = Vector([1.0, 0, 0]) 543 544 self.add_global_vectors("support_normal", "majoraxis") 545 self.add_global_points("centre")
Construct a toroidal mirror.
Parameters
MajorRadius : float
Major radius of the toroid in mm.
MinorRadius : float
Minor radius of the toroid in mm.
Support : ART.ModuleSupport.Support
def
get_local_normal(self, Point):
581 def get_local_normal(self, Point): 582 """Return the normal unit vector on the toroidal surface at point Point in the reference frame of the mirror""" 583 x = Point[0] 584 y = Point[1] 585 z = Point[2] 586 A = self.majorradius**2 - self.minorradius**2 587 588 Gradient = Vector(np.zeros(3)) 589 Gradient[0] = ( 590 4 * (x**3 + x * y**2 + x * z**2 + x * A) 591 - 8 * x * self.majorradius**2 592 ) 593 Gradient[1] = 4 * (y**3 + y * x**2 + y * z**2 + y * A) 594 Gradient[2] = ( 595 4 * (z**3 + z * x**2 + z * y**2 + z * A) 596 - 8 * z * self.majorradius**2 597 ) 598 599 return -Gradient.normalized()
Return the normal unit vector on the toroidal surface at point Point in the reference frame of the mirror
def
render( Mirror, Npoints=1000, draw_support=False, draw_points=False, draw_vectors=True, recenter_support=True):
234def _RenderMirror(Mirror, Npoints=1000, draw_support = False, draw_points = False, draw_vectors = True , recenter_support = True): 235 """ 236 This function renders a mirror in 3D. 237 """ 238 mesh = Mirror._Render(Npoints) 239 p = pv.Plotter() 240 p.add_mesh(mesh) 241 if draw_support: 242 support = Mirror.support._Render() 243 if recenter_support: 244 support.translate(Mirror.r0,inplace=True) 245 p.add_mesh(support, color="gray", opacity=0.5) 246 247 if draw_vectors: 248 # We draw the important vectors of the optical element 249 # For that, if we have a "vectors" attribute, we use that 250 # (a dictionary with the vector names as keys and the colors as values) 251 # Otherwise we use the default vectors: "support_normal", "majoraxis" 252 if hasattr(Mirror, "vectors"): 253 vectors = Mirror.vectors 254 else: 255 vectors = {"support_normal_ref": "red", "majoraxis_ref": "blue"} 256 for vector, color in vectors.items(): 257 if hasattr(Mirror, vector): 258 p.add_arrows(Mirror.r0, 10*getattr(Mirror, vector), color=color) 259 260 if draw_points: 261 # We draw the important points of the optical element 262 # For that, if we have a "points" attribute, we use that 263 # (a dictionary with the point names as keys and the colors as values) 264 # Otherwise we use the default points: "centre_ref" 265 if hasattr(Mirror, "points"): 266 points = Mirror.points 267 else: 268 points = {"centre_ref": "red"} 269 for point, color in points.items(): 270 if hasattr(Mirror, point): 271 p.add_mesh(pv.Sphere(radius=1, center=getattr(Mirror, point)), color=color) 272 else: 273 p.add_mesh(pv.Sphere(radius=1, center=Mirror.r0), color="red") 274 p.show() 275 return p
This function renders a mirror in 3D.
Inherited Members
def
ReturnOptimalToroidalRadii( Focal: float, AngleIncidence: float) -> (<class 'float'>, <class 'float'>):
608def ReturnOptimalToroidalRadii( 609 Focal: float, AngleIncidence: float 610) -> (float, float): 611 """ 612 Get optimal parameters for a toroidal mirror. 613 614 Useful helper function to get the optimal major and minor radii for a toroidal mirror to achieve a 615 focal length 'Focal' with an angle of incidence 'AngleIncidence' and with vanishing astigmatism. 616 617 Parameters 618 ---------- 619 Focal : float 620 Focal length in mm. 621 622 AngleIncidence : int 623 Angle of incidence in degrees. 624 625 Returns 626 ------- 627 OptimalMajorRadius, OptimalMinorRadius : float, float. 628 """ 629 AngleIncidenceRadian = AngleIncidence * np.pi / 180 630 OptimalMajorRadius = ( 631 2 632 * Focal 633 * (1 / np.cos(AngleIncidenceRadian) - np.cos(AngleIncidenceRadian)) 634 ) 635 OptimalMinorRadius = 2 * Focal * np.cos(AngleIncidenceRadian) 636 return OptimalMajorRadius, OptimalMinorRadius
Get optimal parameters for a toroidal mirror.
Useful helper function to get the optimal major and minor radii for a toroidal mirror to achieve a focal length 'Focal' with an angle of incidence 'AngleIncidence' and with vanishing astigmatism.
Parameters
Focal : float
Focal length in mm.
AngleIncidence : int
Angle of incidence in degrees.
Returns
OptimalMajorRadius, OptimalMinorRadius : float, float.
640class MirrorEllipsoidal(Mirror): 641 """ 642 Ellipdoidal mirror surface with eqn. $(x/a)^2 + (y/b)^2 + (z/b)^2 = 1$, where $a$ and $b$ are semi major and semi minor axes. 643 644  645 646 Attributes 647 ---------- 648 a : float 649 Semi major axis of the ellipsoid in mm. 650 651 b : float 652 Semi minor axis of the ellipsoid in mm. 653 654 support : ART.ModuleSupport.Support 655 656 type : str 'Ellipsoidal Mirror'. 657 658 Methods 659 ------- 660 MirrorEllipsoidal.get_normal(Point) 661 662 MirrorEllipsoidal.get_centre() 663 664 MirrorEllipsoidal.get_grid3D(NbPoints) 665 666 """ 667 668 def __init__( 669 self, 670 Support, 671 SemiMajorAxis=None, 672 SemiMinorAxis=None, 673 OffAxisAngle=None, 674 f_object=None, 675 f_image=None, 676 IncidenceAngle=None, 677 Magnification=None, 678 DistanceObjectImage=None, 679 Eccentricity=None, 680 **kwargs, 681 ): 682 """ 683 Generate an ellipsoidal mirror with given parameters. 684 685 The angles are given in degrees but converted to radians for internal calculations. 686 You can specify the semi-major and semi-minor axes, the off-axis angle, the object and image focal distances or some other parameters. 687 The constructor uses a magic calculator to determine the missing parameters. 688 689 Parameters 690 ---------- 691 Support : TYPE 692 ART.ModuleSupport.Support. 693 SemiMajorAxis : float (optional) 694 Semi major axis of the ellipsoid in mm.. 695 SemiMinorAxis : float (optional) 696 Semi minor axis of the ellipsoid in mm.. 697 OffAxisAngle : float (optional) 698 Off-axis angle of the mirror in mm. Defined as the angle at the centre of the mirror between the two foci.. 699 f_object : float (optional) 700 Object focal distance in mm. 701 f_image : float (optional) 702 Image focal distance in mm. 703 IncidenceAngle : float (optional) 704 Angle of incidence in degrees. 705 Magnification : float (optional) 706 Magnification. 707 DistanceObjectImage : float (optional) 708 Distance between object and image in mm. 709 Eccentricity : float (optional) 710 Eccentricity of the ellipsoid. 711 """ 712 super().__init__() 713 self.type = "Ellipsoidal Mirror" 714 self.support = Support 715 716 if "Surface" in kwargs: 717 self.Surface = kwargs["Surface"] 718 719 parameter_calculator = EllipsoidalMirrorCalculator() 720 values, steps = parameter_calculator.calculate_values(verify_consistency=True, 721 f1=f_object, 722 f2=f_image, 723 a=SemiMajorAxis, 724 b=SemiMinorAxis, 725 offset_angle=OffAxisAngle, 726 incidence_angle=IncidenceAngle, 727 m=Magnification, 728 l=DistanceObjectImage, 729 e=Eccentricity) 730 self.a = values["a"] 731 self.b = values["b"] 732 self._offaxisangle = np.deg2rad(values["offset_angle"]) 733 self._f_object = values["f1"] 734 self._f_image = values["f2"] 735 736 self.centre_ref = self._get_centre_ref() 737 738 self.r0 = self.centre_ref 739 740 self.support_normal_ref = Vector([0, 0, 1.0]) 741 self.majoraxis_ref = Vector([1.0, 0, 0]) 742 743 self.f1_ref = Point([-self.a, 0,0]) 744 self.f2_ref = Point([ self.a, 0,0]) 745 self.towards_image_ref = (self.f2_ref - self.centre_ref).normalized() 746 self.towards_object_ref = (self.f1_ref - self.centre_ref).normalized() 747 self.centre_normal_ref = self.get_local_normal(self.centre_ref) 748 749 self.add_global_points("f1", "f2", "centre") 750 self.add_global_vectors("towards_image", "towards_object", "centre_normal", "support_normal", "majoraxis") 751 752 def _get_intersection(self, Ray): 753 """Return the intersection point between Ray and the ellipsoidal mirror surface.""" 754 ux, uy, uz = Ray.vector 755 xA, yA, zA = Ray.point 756 757 da = (uy**2 + uz**2) / self.b**2 + (ux / self.a) ** 2 758 db = 2 * ((uy * yA + uz * zA) / self.b**2 + (ux * xA) / self.a**2) 759 dc = (yA**2 + zA**2) / self.b**2 + (xA / self.a) ** 2 - 1 760 761 Solution = mgeo.SolverQuadratic(da, db, dc) 762 Solution = mgeo.KeepPositiveSolution(Solution) 763 764 ListPointIntersection = [] 765 C = self.get_centre_ref() 766 for t in Solution: 767 Intersect = Ray.vector * t + Ray.point 768 if Intersect[2] < 0 and Intersect - C in self.support: 769 ListPointIntersection.append(Intersect) 770 771 return _IntersectionRayMirror(Ray.point, ListPointIntersection) 772 773 def get_local_normal(self, Point): 774 """Return the normal unit vector on the ellipsoidal surface at point Point.""" 775 Gradient = Vector(np.zeros(3)) 776 777 Gradient[0] = -Point[0] / self.a**2 778 Gradient[1] = -Point[1] / self.b**2 779 Gradient[2] = -Point[2] / self.b**2 780 781 return Gradient.normalized() 782 783 def _get_centre_ref(self): 784 """Return 3D coordinates of the point on the mirror surface at the center of its support.""" 785 foci = 2 * np.sqrt(self.a**2 - self.b**2) # distance between the foci 786 h = -foci / 2 / np.tan(self._offaxisangle) 787 R = np.sqrt(foci**2 / 4 + h**2) 788 sign = 1 789 if math.isclose(self._offaxisangle, np.pi / 2): 790 h = 0 791 elif self._offaxisangle > np.pi / 2: 792 h = -h 793 sign = -1 794 a = 1 - self.a**2 / self.b**2 795 b = -2 * h 796 c = self.a**2 + h**2 - R**2 797 z = (-b + sign * np.sqrt(b**2 - 4 * a * c)) / (2 * a) 798 if math.isclose(z**2, self.b**2): 799 return np.array([0, 0, -self.b]) 800 x = self.a * np.sqrt(1 - z**2 / self.b**2) 801 centre = Point([x, 0, sign * z]) 802 return centre 803 804 def _zfunc(self, PointArray): 805 x = PointArray[:,0] 806 y = PointArray[:,1] 807 x-= self.r0[0] 808 y-= self.r0[1] 809 z = -np.sqrt(1 - (x / self.a)**2 - (y / self.b)**2) 810 z += self.r0[2] 811 return z
Ellipdoidal mirror surface with eqn. $(x/a)^2 + (y/b)^2 + (z/b)^2 = 1$, where $a$ and $b$ are semi major and semi minor axes.
Attributes
a : float
Semi major axis of the ellipsoid in mm.
b : float
Semi minor axis of the ellipsoid in mm.
support : ART.ModuleSupport.Support
type : str 'Ellipsoidal Mirror'.
Methods
MirrorEllipsoidal.get_normal(Point)
MirrorEllipsoidal.get_centre()
MirrorEllipsoidal.get_grid3D(NbPoints)
MirrorEllipsoidal( Support, SemiMajorAxis=None, SemiMinorAxis=None, OffAxisAngle=None, f_object=None, f_image=None, IncidenceAngle=None, Magnification=None, DistanceObjectImage=None, Eccentricity=None, **kwargs)
668 def __init__( 669 self, 670 Support, 671 SemiMajorAxis=None, 672 SemiMinorAxis=None, 673 OffAxisAngle=None, 674 f_object=None, 675 f_image=None, 676 IncidenceAngle=None, 677 Magnification=None, 678 DistanceObjectImage=None, 679 Eccentricity=None, 680 **kwargs, 681 ): 682 """ 683 Generate an ellipsoidal mirror with given parameters. 684 685 The angles are given in degrees but converted to radians for internal calculations. 686 You can specify the semi-major and semi-minor axes, the off-axis angle, the object and image focal distances or some other parameters. 687 The constructor uses a magic calculator to determine the missing parameters. 688 689 Parameters 690 ---------- 691 Support : TYPE 692 ART.ModuleSupport.Support. 693 SemiMajorAxis : float (optional) 694 Semi major axis of the ellipsoid in mm.. 695 SemiMinorAxis : float (optional) 696 Semi minor axis of the ellipsoid in mm.. 697 OffAxisAngle : float (optional) 698 Off-axis angle of the mirror in mm. Defined as the angle at the centre of the mirror between the two foci.. 699 f_object : float (optional) 700 Object focal distance in mm. 701 f_image : float (optional) 702 Image focal distance in mm. 703 IncidenceAngle : float (optional) 704 Angle of incidence in degrees. 705 Magnification : float (optional) 706 Magnification. 707 DistanceObjectImage : float (optional) 708 Distance between object and image in mm. 709 Eccentricity : float (optional) 710 Eccentricity of the ellipsoid. 711 """ 712 super().__init__() 713 self.type = "Ellipsoidal Mirror" 714 self.support = Support 715 716 if "Surface" in kwargs: 717 self.Surface = kwargs["Surface"] 718 719 parameter_calculator = EllipsoidalMirrorCalculator() 720 values, steps = parameter_calculator.calculate_values(verify_consistency=True, 721 f1=f_object, 722 f2=f_image, 723 a=SemiMajorAxis, 724 b=SemiMinorAxis, 725 offset_angle=OffAxisAngle, 726 incidence_angle=IncidenceAngle, 727 m=Magnification, 728 l=DistanceObjectImage, 729 e=Eccentricity) 730 self.a = values["a"] 731 self.b = values["b"] 732 self._offaxisangle = np.deg2rad(values["offset_angle"]) 733 self._f_object = values["f1"] 734 self._f_image = values["f2"] 735 736 self.centre_ref = self._get_centre_ref() 737 738 self.r0 = self.centre_ref 739 740 self.support_normal_ref = Vector([0, 0, 1.0]) 741 self.majoraxis_ref = Vector([1.0, 0, 0]) 742 743 self.f1_ref = Point([-self.a, 0,0]) 744 self.f2_ref = Point([ self.a, 0,0]) 745 self.towards_image_ref = (self.f2_ref - self.centre_ref).normalized() 746 self.towards_object_ref = (self.f1_ref - self.centre_ref).normalized() 747 self.centre_normal_ref = self.get_local_normal(self.centre_ref) 748 749 self.add_global_points("f1", "f2", "centre") 750 self.add_global_vectors("towards_image", "towards_object", "centre_normal", "support_normal", "majoraxis")
Generate an ellipsoidal mirror with given parameters.
The angles are given in degrees but converted to radians for internal calculations. You can specify the semi-major and semi-minor axes, the off-axis angle, the object and image focal distances or some other parameters. The constructor uses a magic calculator to determine the missing parameters.
Parameters
Support : TYPE ART.ModuleSupport.Support. SemiMajorAxis : float (optional) Semi major axis of the ellipsoid in mm.. SemiMinorAxis : float (optional) Semi minor axis of the ellipsoid in mm.. OffAxisAngle : float (optional) Off-axis angle of the mirror in mm. Defined as the angle at the centre of the mirror between the two foci.. f_object : float (optional) Object focal distance in mm. f_image : float (optional) Image focal distance in mm. IncidenceAngle : float (optional) Angle of incidence in degrees. Magnification : float (optional) Magnification. DistanceObjectImage : float (optional) Distance between object and image in mm. Eccentricity : float (optional) Eccentricity of the ellipsoid.
def
get_local_normal(self, Point):
773 def get_local_normal(self, Point): 774 """Return the normal unit vector on the ellipsoidal surface at point Point.""" 775 Gradient = Vector(np.zeros(3)) 776 777 Gradient[0] = -Point[0] / self.a**2 778 Gradient[1] = -Point[1] / self.b**2 779 Gradient[2] = -Point[2] / self.b**2 780 781 return Gradient.normalized()
Return the normal unit vector on the ellipsoidal surface at point Point.
def
render( Mirror, Npoints=1000, draw_support=False, draw_points=False, draw_vectors=True, recenter_support=True):
234def _RenderMirror(Mirror, Npoints=1000, draw_support = False, draw_points = False, draw_vectors = True , recenter_support = True): 235 """ 236 This function renders a mirror in 3D. 237 """ 238 mesh = Mirror._Render(Npoints) 239 p = pv.Plotter() 240 p.add_mesh(mesh) 241 if draw_support: 242 support = Mirror.support._Render() 243 if recenter_support: 244 support.translate(Mirror.r0,inplace=True) 245 p.add_mesh(support, color="gray", opacity=0.5) 246 247 if draw_vectors: 248 # We draw the important vectors of the optical element 249 # For that, if we have a "vectors" attribute, we use that 250 # (a dictionary with the vector names as keys and the colors as values) 251 # Otherwise we use the default vectors: "support_normal", "majoraxis" 252 if hasattr(Mirror, "vectors"): 253 vectors = Mirror.vectors 254 else: 255 vectors = {"support_normal_ref": "red", "majoraxis_ref": "blue"} 256 for vector, color in vectors.items(): 257 if hasattr(Mirror, vector): 258 p.add_arrows(Mirror.r0, 10*getattr(Mirror, vector), color=color) 259 260 if draw_points: 261 # We draw the important points of the optical element 262 # For that, if we have a "points" attribute, we use that 263 # (a dictionary with the point names as keys and the colors as values) 264 # Otherwise we use the default points: "centre_ref" 265 if hasattr(Mirror, "points"): 266 points = Mirror.points 267 else: 268 points = {"centre_ref": "red"} 269 for point, color in points.items(): 270 if hasattr(Mirror, point): 271 p.add_mesh(pv.Sphere(radius=1, center=getattr(Mirror, point)), color=color) 272 else: 273 p.add_mesh(pv.Sphere(radius=1, center=Mirror.r0), color="red") 274 p.show() 275 return p
This function renders a mirror in 3D.
vectors =
{'support_normal_ref': 'red', 'majoraxis_ref': 'blue', 'towards_image_ref': 'green', 'towards_object_ref': 'green', 'centre_normal_ref': 'purple'}
Inherited Members
815class MirrorCylindrical(Mirror): 816 """ 817 Cylindrical mirror surface with eqn. $y^2 + z^2 = R^2$, where $R$ is the radius. 818 819 Attributes 820 ---------- 821 radius : float 822 Radius of curvature. A postive value for concave mirror, a negative value for a convex mirror. 823 824 support : ART.ModuleSupport.Support 825 826 type : str 'Cylindrical Mirror'. 827 828 Methods 829 ------- 830 MirrorCylindrical.get_normal(Point) 831 832 MirrorCylindrical.get_centre() 833 834 MirrorCylindrical.get_grid3D(NbPoints) 835 """ 836 837 def __init__(self, Support, Radius, **kwargs): 838 """ 839 Construct a cylindrical mirror. 840 841 Parameters 842 ---------- 843 Radius : float 844 The radius of curvature in mm. A postive value for concave mirror, a negative value for a convex mirror. 845 846 Support : ART.ModuleSupport.Support 847 848 """ 849 super().__init__() 850 if Radius < 0: 851 self.type = "CylindricalCX Mirror" 852 self.curvature = Curvature.CONVEX 853 self.radius = -Radius 854 else: 855 self.type = "CylindricalCC Mirror" 856 self.curvature = Curvature.CONCAVE 857 self.radius = Radius 858 859 self.support = Support 860 861 if "Surface" in kwargs: 862 self.Surface = kwargs["Surface"] 863 864 self.r0 = Point([0.0, 0.0, -self.radius]) 865 866 self.support_normal_ref = Vector([0, 0, 1.0]) 867 self.majoraxis_ref = Vector([1.0, 0, 0]) 868 self.centre_ref = Point([0.0, 0.0, -self.radius]) 869 870 self.add_global_points("centre") 871 self.add_global_vectors("support_normal", "majoraxis") 872 873 def _get_intersection(self, Ray): 874 """Return the intersection point between the Ray and the cylinder.""" 875 uy = Ray.vector[1] 876 uz = Ray.vector[2] 877 yA = Ray.point[1] 878 zA = Ray.point[2] 879 880 a = uy**2 + uz**2 881 b = 2 * (uy * yA + uz * zA) 882 c = yA**2 + zA**2 - self.radius**2 883 884 Solution = mgeo.SolverQuadratic(a, b, c) 885 Solution = mgeo.KeepPositiveSolution(Solution) 886 887 ListPointIntersection = [] 888 for t in Solution: 889 Intersect = Ray.vector * t + Ray.point 890 if Intersect[2] < 0 and Intersect in self.support: 891 ListPointIntersection.append(Intersect) 892 893 return _IntersectionRayMirror(Ray.point, ListPointIntersection) 894 895 def get_local_normal(self, Point): 896 """Return the normal unit vector on the cylinder surface at point P.""" 897 Gradient = Vector([0, -Point[1], -Point[2]]) 898 return Gradient.normalized() 899 900 def get_centre(self): 901 """Return 3D coordinates of the point on the mirror surface at the center of its support.""" 902 return Point([0, 0, -self.radius]) 903 904 def _zfunc(self, PointArray): 905 y = PointArray[:,1] 906 return -np.sqrt(self.radius**2 - y**2)
Cylindrical mirror surface with eqn. $y^2 + z^2 = R^2$, where $R$ is the radius.
Attributes
radius : float
Radius of curvature. A postive value for concave mirror, a negative value for a convex mirror.
support : ART.ModuleSupport.Support
type : str 'Cylindrical Mirror'.
Methods
MirrorCylindrical.get_normal(Point)
MirrorCylindrical.get_centre()
MirrorCylindrical.get_grid3D(NbPoints)
MirrorCylindrical(Support, Radius, **kwargs)
837 def __init__(self, Support, Radius, **kwargs): 838 """ 839 Construct a cylindrical mirror. 840 841 Parameters 842 ---------- 843 Radius : float 844 The radius of curvature in mm. A postive value for concave mirror, a negative value for a convex mirror. 845 846 Support : ART.ModuleSupport.Support 847 848 """ 849 super().__init__() 850 if Radius < 0: 851 self.type = "CylindricalCX Mirror" 852 self.curvature = Curvature.CONVEX 853 self.radius = -Radius 854 else: 855 self.type = "CylindricalCC Mirror" 856 self.curvature = Curvature.CONCAVE 857 self.radius = Radius 858 859 self.support = Support 860 861 if "Surface" in kwargs: 862 self.Surface = kwargs["Surface"] 863 864 self.r0 = Point([0.0, 0.0, -self.radius]) 865 866 self.support_normal_ref = Vector([0, 0, 1.0]) 867 self.majoraxis_ref = Vector([1.0, 0, 0]) 868 self.centre_ref = Point([0.0, 0.0, -self.radius]) 869 870 self.add_global_points("centre") 871 self.add_global_vectors("support_normal", "majoraxis")
Construct a cylindrical mirror.
Parameters
Radius : float
The radius of curvature in mm. A postive value for concave mirror, a negative value for a convex mirror.
Support : ART.ModuleSupport.Support
def
get_local_normal(self, Point):
895 def get_local_normal(self, Point): 896 """Return the normal unit vector on the cylinder surface at point P.""" 897 Gradient = Vector([0, -Point[1], -Point[2]]) 898 return Gradient.normalized()
Return the normal unit vector on the cylinder surface at point P.
def
get_centre(self):
900 def get_centre(self): 901 """Return 3D coordinates of the point on the mirror surface at the center of its support.""" 902 return Point([0, 0, -self.radius])
Return 3D coordinates of the point on the mirror surface at the center of its support.
def
render( Mirror, Npoints=1000, draw_support=False, draw_points=False, draw_vectors=True, recenter_support=True):
234def _RenderMirror(Mirror, Npoints=1000, draw_support = False, draw_points = False, draw_vectors = True , recenter_support = True): 235 """ 236 This function renders a mirror in 3D. 237 """ 238 mesh = Mirror._Render(Npoints) 239 p = pv.Plotter() 240 p.add_mesh(mesh) 241 if draw_support: 242 support = Mirror.support._Render() 243 if recenter_support: 244 support.translate(Mirror.r0,inplace=True) 245 p.add_mesh(support, color="gray", opacity=0.5) 246 247 if draw_vectors: 248 # We draw the important vectors of the optical element 249 # For that, if we have a "vectors" attribute, we use that 250 # (a dictionary with the vector names as keys and the colors as values) 251 # Otherwise we use the default vectors: "support_normal", "majoraxis" 252 if hasattr(Mirror, "vectors"): 253 vectors = Mirror.vectors 254 else: 255 vectors = {"support_normal_ref": "red", "majoraxis_ref": "blue"} 256 for vector, color in vectors.items(): 257 if hasattr(Mirror, vector): 258 p.add_arrows(Mirror.r0, 10*getattr(Mirror, vector), color=color) 259 260 if draw_points: 261 # We draw the important points of the optical element 262 # For that, if we have a "points" attribute, we use that 263 # (a dictionary with the point names as keys and the colors as values) 264 # Otherwise we use the default points: "centre_ref" 265 if hasattr(Mirror, "points"): 266 points = Mirror.points 267 else: 268 points = {"centre_ref": "red"} 269 for point, color in points.items(): 270 if hasattr(Mirror, point): 271 p.add_mesh(pv.Sphere(radius=1, center=getattr(Mirror, point)), color=color) 272 else: 273 p.add_mesh(pv.Sphere(radius=1, center=Mirror.r0), color="red") 274 p.show() 275 return p
This function renders a mirror in 3D.
Inherited Members
913class GrazingParabola(Mirror): 914 r""" 915 A parabolic mirror with a support parallel to the surface at the optical center. 916 Needs to be completed, currently not functional. 917 TODO 918 """ 919 920 def __init__(self, Support, 921 FocalEffective: float=None, 922 OffAxisAngle: float = None, 923 FocalParent: float = None, 924 RadiusParent: float = None, 925 OffAxisDistance: float = None, 926 MoreThan90: bool = None, 927 **kwargs): 928 """ 929 Initialise a Parabolic mirror. 930 931 Parameters 932 ---------- 933 FocalEffective : float 934 Effective focal length of the parabola in mm. 935 936 OffAxisAngle : float 937 Off-axis angle *in degrees* of the parabola. 938 939 Support : ART.ModuleSupport.Support 940 941 """ 942 super().__init__() 943 self.curvature = Curvature.CONCAVE 944 self.support = Support 945 self.type = "Grazing Parabolic Mirror" 946 947 if "Surface" in kwargs: 948 self.Surface = kwargs["Surface"] 949 950 parameter_calculator = OAP_calculator() 951 values,steps = parameter_calculator.calculate_values(verify_consistency=True, 952 fs=FocalEffective, 953 theta=OffAxisAngle, 954 fp=FocalParent, 955 Rc=RadiusParent, 956 OAD=OffAxisDistance, 957 more_than_90=MoreThan90) 958 self._offaxisangle = values["theta"] 959 self._feff = values["fs"] 960 self._p = values["p"] 961 self._offaxisdistance = values["OAD"] 962 self._fparent = values["fp"] 963 self._rparent = values["Rc"] 964 965 #self.r0 = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 966 self.r0 = Point([0, 0, 0]) 967 968 self.centre_ref = Point([0, 0, 0]) 969 self.support_normal_ref = Vector([0, 0, 1.0]) 970 self.majoraxis_ref = Vector([1.0, 0, 0]) 971 972 focus_ref = Point([np.sqrt(self._feff**2 - self._offaxisdistance**2), 0, self._offaxisdistance]) # frame where x is as collimated direction 973 # We need to rotate it in the trigonometric direction by the 90°-theta/2 974 angle = np.pi/2 - self._offaxisangle/2 975 self.focus_ref = Point([focus_ref[0]*np.cos(angle)+focus_ref[2]*np.sin(angle),focus_ref[1], -focus_ref[0]*np.sin(angle)+focus_ref[2]*np.cos(angle), ]) 976 977 self.towards_focusing_ref = (self.focus_ref-self.centre_ref).normalized() 978 towards_collimated_ref = Vector([-1.0, 0, 0]) # Again, we need to rotate it by the same angle 979 self.collimated_ref = Point([-np.cos(angle), 0, np.sin(angle)]) 980 981 self.add_global_points("focus", "centre") 982 self.add_global_vectors("towards_focusing", "towards_collimated", "support_normal", "majoraxis") 983 984 def _get_intersection(self, Ray): 985 """ 986 Return the intersection point between the ray and the parabola. 987 """ 988 pass 989 990 def get_local_normal(self, Point): 991 """ 992 Return the normal unit vector on the paraboloid surface at point Point. 993 """ 994 pass 995 996 def _zfunc(self, PointArray): 997 pass
A parabolic mirror with a support parallel to the surface at the optical center. Needs to be completed, currently not functional. TODO
GrazingParabola( Support, FocalEffective: float = None, OffAxisAngle: float = None, FocalParent: float = None, RadiusParent: float = None, OffAxisDistance: float = None, MoreThan90: bool = None, **kwargs)
920 def __init__(self, Support, 921 FocalEffective: float=None, 922 OffAxisAngle: float = None, 923 FocalParent: float = None, 924 RadiusParent: float = None, 925 OffAxisDistance: float = None, 926 MoreThan90: bool = None, 927 **kwargs): 928 """ 929 Initialise a Parabolic mirror. 930 931 Parameters 932 ---------- 933 FocalEffective : float 934 Effective focal length of the parabola in mm. 935 936 OffAxisAngle : float 937 Off-axis angle *in degrees* of the parabola. 938 939 Support : ART.ModuleSupport.Support 940 941 """ 942 super().__init__() 943 self.curvature = Curvature.CONCAVE 944 self.support = Support 945 self.type = "Grazing Parabolic Mirror" 946 947 if "Surface" in kwargs: 948 self.Surface = kwargs["Surface"] 949 950 parameter_calculator = OAP_calculator() 951 values,steps = parameter_calculator.calculate_values(verify_consistency=True, 952 fs=FocalEffective, 953 theta=OffAxisAngle, 954 fp=FocalParent, 955 Rc=RadiusParent, 956 OAD=OffAxisDistance, 957 more_than_90=MoreThan90) 958 self._offaxisangle = values["theta"] 959 self._feff = values["fs"] 960 self._p = values["p"] 961 self._offaxisdistance = values["OAD"] 962 self._fparent = values["fp"] 963 self._rparent = values["Rc"] 964 965 #self.r0 = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 966 self.r0 = Point([0, 0, 0]) 967 968 self.centre_ref = Point([0, 0, 0]) 969 self.support_normal_ref = Vector([0, 0, 1.0]) 970 self.majoraxis_ref = Vector([1.0, 0, 0]) 971 972 focus_ref = Point([np.sqrt(self._feff**2 - self._offaxisdistance**2), 0, self._offaxisdistance]) # frame where x is as collimated direction 973 # We need to rotate it in the trigonometric direction by the 90°-theta/2 974 angle = np.pi/2 - self._offaxisangle/2 975 self.focus_ref = Point([focus_ref[0]*np.cos(angle)+focus_ref[2]*np.sin(angle),focus_ref[1], -focus_ref[0]*np.sin(angle)+focus_ref[2]*np.cos(angle), ]) 976 977 self.towards_focusing_ref = (self.focus_ref-self.centre_ref).normalized() 978 towards_collimated_ref = Vector([-1.0, 0, 0]) # Again, we need to rotate it by the same angle 979 self.collimated_ref = Point([-np.cos(angle), 0, np.sin(angle)]) 980 981 self.add_global_points("focus", "centre") 982 self.add_global_vectors("towards_focusing", "towards_collimated", "support_normal", "majoraxis")
Initialise a Parabolic mirror.
Parameters
FocalEffective : float
Effective focal length of the parabola in mm.
OffAxisAngle : float
Off-axis angle *in degrees* of the parabola.
Support : ART.ModuleSupport.Support
def
get_local_normal(self, Point):
990 def get_local_normal(self, Point): 991 """ 992 Return the normal unit vector on the paraboloid surface at point Point. 993 """ 994 pass
Return the normal unit vector on the paraboloid surface at point Point.
Inherited Members
def
ReflectionMirrorRayList(Mirror, ListRay, IgnoreDefects=False):
1033def ReflectionMirrorRayList(Mirror, ListRay, IgnoreDefects=False): 1034 """ 1035 Return the the reflected rays according to the law of reflection for the list of incident rays ListRay. 1036 1037 Rays that do not hit the support are not further propagated. 1038 1039 Updates the reflected rays' incidence angle and path. 1040 1041 Parameters 1042 ---------- 1043 Mirror : Mirror-object 1044 1045 ListRay : list[Ray-object] 1046 1047 """ 1048 Deformed = type(Mirror) == DeformedMirror 1049 ListRayReflected = [] 1050 for k in ListRay: 1051 PointMirror = Mirror._get_intersection(k) 1052 1053 if PointMirror is not None: 1054 if Deformed and IgnoreDefects: 1055 M = Mirror.Mirror 1056 else: 1057 M = Mirror 1058 RayReflected = _ReflectionMirrorRay(M, PointMirror, k) 1059 ListRayReflected.append(RayReflected) 1060 return ListRayReflected
Return the the reflected rays according to the law of reflection for the list of incident rays ListRay.
Rays that do not hit the support are not further propagated.
Updates the reflected rays' incidence angle and path.
Parameters
Mirror : Mirror-object
ListRay : list[Ray-object]
1064class DeformedMirror(Mirror): 1065 def __init__(self, Mirror, DeformationList): 1066 self.Mirror = Mirror 1067 self.DeformationList = DeformationList 1068 self.type = Mirror.type 1069 self.support = self.Mirror.support 1070 1071 def get_local_normal(self, PointMirror): 1072 base_normal = self.Mirror.get_normal(PointMirror) 1073 C = self.get_centre() 1074 defects_normals = [ 1075 d.get_normal(PointMirror - C) for d in self.DeformationList 1076 ] 1077 for i in defects_normals: 1078 base_normal = mgeo.normal_add(base_normal, i) 1079 base_normal /= np.linalg.norm(base_normal) 1080 return base_normal 1081 1082 def get_centre(self): 1083 return self.Mirror.get_centre() 1084 1085 def get_grid3D(self, NbPoint, **kwargs): 1086 return self.Mirror.get_grid3D(NbPoint, **kwargs) 1087 1088 def _get_intersection(self, Ray): 1089 Intersect = self.Mirror._get_intersection(Ray) 1090 if Intersect is not None: 1091 h = sum( 1092 D.get_offset(Intersect - self.get_centre()) 1093 for D in self.DeformationList 1094 ) 1095 alpha = mgeo.AngleBetweenTwoVectors( 1096 -Ray.vector, self.Mirror.get_normal(Intersect) 1097 ) 1098 Intersect -= Ray.vector * h / np.cos(alpha) 1099 return Intersect
Abstract base class for mirrors.
def
get_local_normal(self, PointMirror):
1071 def get_local_normal(self, PointMirror): 1072 base_normal = self.Mirror.get_normal(PointMirror) 1073 C = self.get_centre() 1074 defects_normals = [ 1075 d.get_normal(PointMirror - C) for d in self.DeformationList 1076 ] 1077 for i in defects_normals: 1078 base_normal = mgeo.normal_add(base_normal, i) 1079 base_normal /= np.linalg.norm(base_normal) 1080 return base_normal
This method should return the normal unit vector in point 'Point' on the mirror surface. The normal is in the reference frame of the mirror.