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 for t in Solution: 288 Intersect = Ray.vector * t + Ray.point 289 if Intersect - self.r0 in self.support: 290 return Intersect, True 291 292 return Ray.point, False 293 294 def get_local_normal(self, Point): 295 """Return the normal unit vector on the spherical surface at point Point.""" 296 return -Vector(Point).normalized() 297 298 def _zfunc(self, PointArray): 299 x = PointArray[:,0] 300 y = PointArray[:,1] 301 return -np.sqrt(self.radius**2 - x**2 - y**2) 302 303# %% Parabolic mirror definitions 304class MirrorParabolic(Mirror): 305 r""" 306 A paraboloid with vertex at the origin $O=[0,0,0]$ and symmetry axis z: 307 $z = \frac{1}{4f}[x^2 + y^2]$ where $f$ is the focal lenght of the *mother* 308 parabola (i.e. measured from its center at $O$ to the focal point $F$). 309 310 The center of the support is shifted along the x-direction by the off-axis distance $x_c$. 311 This leads to an *effective focal length* $f_\mathrm{eff}$, measured from the shifted center 312 of the support $P$ to the focal point $F$. 313 It is related to the mother focal length by $f = f_\\mathrm{eff} \cos^2(\alpha/2) $, 314 or equivalently $ p = 2f = f_\mathrm{eff} (1+\\cos\\alpha)$, where $\\alpha$ 315 is the off-axis angle, and $p = 2f$ is called the semi latus rectum. 316 317 Another useful relationship is that between the off-axis distance and the resulting 318 off-axis angle: $x_c = 2 f \tan(\alpha/2)$. 319 320 321  322 323 Attributes 324 ---------- 325 offaxisangle : float 326 Off-axis angle of the parabola. Modifying this also updates p, keeping feff constant. 327 Attention: The off-axis angle must be *given in degrees*, but is stored and will be *returned in radian* ! 328 329 feff : float 330 Effective focal length of the parabola in mm. Modifying this also updates p, keeping offaxisangle constant. 331 332 p : float 333 Semi latus rectum of the parabola in mm. Modifying this also updates feff, keeping offaxisangle constant. 334 335 support : ART.ModuleSupport.Support 336 337 type : str 'Parabolic Mirror'. 338 339 Methods 340 ------- 341 MirrorParabolic.get_normal(Point) 342 343 MirrorParabolic.get_centre() 344 345 MirrorParabolic.get_grid3D(NbPoints) 346 347 """ 348 #vectorised = True # Unitl I fix the closest solution thing 349 def __init__(self, Support, 350 FocalEffective: float=None, 351 OffAxisAngle: float = None, 352 FocalParent: float = None, 353 RadiusParent: float = None, 354 OffAxisDistance: float = None, 355 MoreThan90: bool = None, 356 **kwargs): 357 """ 358 Initialise a Parabolic mirror. 359 360 Parameters 361 ---------- 362 FocalEffective : float 363 Effective focal length of the parabola in mm. 364 365 OffAxisAngle : float 366 Off-axis angle *in degrees* of the parabola. 367 368 Support : ART.ModuleSupport.Support 369 370 """ 371 super().__init__() 372 self.curvature = Curvature.CONCAVE 373 self.support = Support 374 self.type = "Parabolic Mirror" 375 376 if "Surface" in kwargs: 377 self.Surface = kwargs["Surface"] 378 379 parameter_calculator = OAP_calculator() 380 values,steps = parameter_calculator.calculate_values(verify_consistency=True, 381 fs=FocalEffective, 382 theta=OffAxisAngle, 383 fp=FocalParent, 384 Rc=RadiusParent, 385 OAD=OffAxisDistance, 386 more_than_90=MoreThan90) 387 self._offaxisangle = values["theta"] 388 self._feff = values["fs"] 389 self._p = values["p"] 390 self._offaxisdistance = values["OAD"] 391 self._fparent = values["fp"] 392 self._rparent = values["Rc"] 393 394 self.r0 = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 395 396 self.centre_ref = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 397 self.support_normal_ref = Vector([0, 0, 1.0]) 398 self.majoraxis_ref = Vector([1.0, 0, 0]) 399 400 self.focus_ref = Point([0.0, 0, self._fparent]) 401 self.towards_focusing_ref = (self.focus_ref-self.centre_ref).normalized() 402 self.towards_collimated_ref = Vector([0, 0, 1.0]) 403 404 self.add_global_points("focus", "centre") 405 self.add_global_vectors("towards_focusing", "towards_collimated", "support_normal", "majoraxis") 406 407 def _get_intersection(self, Ray): 408 """Return the intersection point between the ray and the parabola.""" 409 ux, uy, uz = Ray.vector 410 xA, yA, zA = Ray.point 411 412 da = ux**2 + uy**2 413 db = 2 * (ux * xA + uy * yA) - 2 * self._p * uz 414 dc = xA**2 + yA**2 - 2 * self._p * zA 415 416 Solution = mgeo.SolverQuadratic(da, db, dc) 417 Solution = mgeo.KeepPositiveSolution(Solution) 418 IntersectionPoint = [Ray.point + Ray.vector * t for t in Solution if t > 0] 419 distances = [mgeo.Vector(i - self.r0).norm for i in IntersectionPoint] 420 IntersectionPoint = IntersectionPoint[distances.index(min(distances))] 421 return IntersectionPoint, IntersectionPoint-self.r0 in self.support 422 423 def _get_intersections(self, RayList): 424 """ 425 Vectorised version of the intersection calculation. 426 """ 427 vectors = np.array([ray.vector for ray in RayList]) 428 points = np.array([ray.point for ray in RayList]) 429 Solutions = np.zeros(len(RayList), dtype=float) 430 431 da = np.sum(vectors[:,0:2]**2, axis=1) 432 db = 2 * (np.sum(vectors[:,0:2] * points[:,0:2], axis=1) - self._p * vectors[:,2]) 433 dc = np.sum(points[:,0:2]**2, axis=1) - 2 * self._p * points[:,2] 434 435 linear = da == 0 436 Solutions[linear] = -dc[linear] / db[linear] 437 438 centered = (dc==0) & (da!=0) # If equation is of form ax**2 + bx = 0 439 Solutions[centered] = np.maximum(-db[centered] / da[centered], 0) 440 441 unstable = (4*da*dc) < db**2 * 1e-10 # Cases leading to catastrophic cancellation 442 443 Deltas = db**2 - 4 * da * dc 444 OK = Deltas >= 0 445 SolutionsPlus = (-db + np.sqrt(Deltas)) / 2 / da 446 SolutionsMinus = (-db - np.sqrt(Deltas)) / 2 / da 447 SolutionsPlus[unstable] = dc[unstable] / (da[unstable] * SolutionsMinus[unstable]) 448 SolutionsPlus = np.where(SolutionsPlus >= 0, SolutionsPlus, np.inf) 449 SolutionsMinus = np.where(SolutionsMinus >= 0, SolutionsMinus, np.inf) 450 Solutions = np.minimum(SolutionsPlus, SolutionsMinus) 451 Solutions = np.maximum(Solutions, 0) 452 Points = points + vectors * Solutions[:,np.newaxis] 453 OK = OK & (Solutions > 0) & np.array([p-self.r0 in self.support for p in Points]) 454 return mgeo.PointArray(Points), OK 455 456 def get_local_normal(self, Point): 457 """Return the normal unit vector on the paraboloid surface at point Point.""" 458 Gradient = Vector(np.zeros(3)) 459 Gradient[0] = -Point[0] 460 Gradient[1] = -Point[1] 461 Gradient[2] = self._p 462 return Gradient.normalized() 463 464 def get_local_normals(self, Points): 465 """ 466 Vectorised version of the normal calculation. 467 """ 468 Gradients = mgeo.VectorArray(np.zeros_like(Points)) 469 Gradients[:,0] = -Points[:,0] 470 Gradients[:,1] = -Points[:,1] 471 Gradients[:,2] = self._p 472 return Gradients.normalized() 473 474 def _zfunc(self, PointArray): 475 x = PointArray[:,0] 476 y = PointArray[:,1] 477 return (x**2 + y**2) / 2 / self._p 478 479 480# %% Toroidal mirror definitions 481class MirrorToroidal(Mirror): 482 r""" 483 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. 484 485  486 487 Attributes 488 ---------- 489 majorradius : float 490 Major radius of the toroid in mm. 491 492 minorradius : float 493 Minor radius of the toroid in mm. 494 495 support : ART.ModuleSupport.Support 496 497 type : str 'Toroidal Mirror'. 498 499 Methods 500 ------- 501 MirrorToroidal.get_normal(Point) 502 503 MirrorToroidal.get_centre() 504 505 MirrorToroidal.get_grid3D(NbPoints) 506 507 """ 508 509 def __init__(self, Support, MajorRadius, MinorRadius, **kwargs): 510 """ 511 Construct a toroidal mirror. 512 513 Parameters 514 ---------- 515 MajorRadius : float 516 Major radius of the toroid in mm. 517 518 MinorRadius : float 519 Minor radius of the toroid in mm. 520 521 Support : ART.ModuleSupport.Support 522 523 """ 524 super().__init__() 525 self.support = Support 526 self.type = "Toroidal Mirror" 527 528 if "Surface" in kwargs: 529 self.Surface = kwargs["Surface"] 530 531 self.curvature = Curvature.CONCAVE 532 533 self.majorradius = MajorRadius 534 self.minorradius = MinorRadius 535 536 self.r0 = Point([0.0, 0.0, -MajorRadius - MinorRadius]) 537 538 self.centre_ref = Point([0.0, 0.0, -MajorRadius - MinorRadius]) 539 self.support_normal_ref = Vector([0, 0, 1.0]) 540 self.majoraxis_ref = Vector([1.0, 0, 0]) 541 542 self.add_global_vectors("support_normal", "majoraxis") 543 self.add_global_points("centre") 544 545 def _get_intersection(self, Ray): 546 """Return the intersection point between Ray and the toroidal mirror surface. The intersection points are given in the reference frame of the mirror""" 547 ux = Ray.vector[0] 548 uz = Ray.vector[2] 549 xA = Ray.point[0] 550 zA = Ray.point[2] 551 552 G = 4.0 * self.majorradius**2 * (ux**2 + uz**2) 553 H = 8.0 * self.majorradius**2 * (ux * xA + uz * zA) 554 I = 4.0 * self.majorradius**2 * (xA**2 + zA**2) 555 J = np.dot(Ray.vector, Ray.vector) 556 K = 2.0 * np.dot(Ray.vector, Ray.point) 557 L = ( 558 np.dot(Ray.point, Ray.point) 559 + self.majorradius**2 560 - self.minorradius**2 561 ) 562 563 a = J**2 564 b = 2 * J * K 565 c = 2 * J * L + K**2 - G 566 d = 2 * K * L - H 567 e = L**2 - I 568 569 Solution = mgeo.SolverQuartic(a, b, c, d, e) 570 if len(Solution) == 0: 571 return None, False 572 Solution = [t for t in Solution if t > 0] 573 IntersectionPoint = [Ray.point + Ray.vector *i for i in Solution] 574 distances = [mgeo.Vector(i - self.r0).norm for i in IntersectionPoint] 575 IntersectionPoint = IntersectionPoint[distances.index(min(distances))] 576 OK = IntersectionPoint - self.r0 in self.support and np.abs(IntersectionPoint[2]-self.r0[2]) < self.majorradius 577 return IntersectionPoint, OK 578 579 def get_local_normal(self, Point): 580 """Return the normal unit vector on the toroidal surface at point Point in the reference frame of the mirror""" 581 x = Point[0] 582 y = Point[1] 583 z = Point[2] 584 A = self.majorradius**2 - self.minorradius**2 585 586 Gradient = Vector(np.zeros(3)) 587 Gradient[0] = ( 588 4 * (x**3 + x * y**2 + x * z**2 + x * A) 589 - 8 * x * self.majorradius**2 590 ) 591 Gradient[1] = 4 * (y**3 + y * x**2 + y * z**2 + y * A) 592 Gradient[2] = ( 593 4 * (z**3 + z * x**2 + z * y**2 + z * A) 594 - 8 * z * self.majorradius**2 595 ) 596 597 return -Gradient.normalized() 598 599 def _zfunc(self, PointArray): 600 x = PointArray[:,0] 601 y = PointArray[:,1] 602 return -np.sqrt( 603 (np.sqrt(self.minorradius**2 - y**2) + self.majorradius)**2 - x**2 604 ) 605 606def ReturnOptimalToroidalRadii( 607 Focal: float, AngleIncidence: float 608) -> (float, float): 609 """ 610 Get optimal parameters for a toroidal mirror. 611 612 Useful helper function to get the optimal major and minor radii for a toroidal mirror to achieve a 613 focal length 'Focal' with an angle of incidence 'AngleIncidence' and with vanishing astigmatism. 614 615 Parameters 616 ---------- 617 Focal : float 618 Focal length in mm. 619 620 AngleIncidence : int 621 Angle of incidence in degrees. 622 623 Returns 624 ------- 625 OptimalMajorRadius, OptimalMinorRadius : float, float. 626 """ 627 AngleIncidenceRadian = AngleIncidence * np.pi / 180 628 OptimalMajorRadius = ( 629 2 630 * Focal 631 * (1 / np.cos(AngleIncidenceRadian) - np.cos(AngleIncidenceRadian)) 632 ) 633 OptimalMinorRadius = 2 * Focal * np.cos(AngleIncidenceRadian) 634 return OptimalMajorRadius, OptimalMinorRadius 635 636 637# %% Ellipsoidal mirror definitions 638class MirrorEllipsoidal(Mirror): 639 """ 640 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. 641 642  643 644 Attributes 645 ---------- 646 a : float 647 Semi major axis of the ellipsoid in mm. 648 649 b : float 650 Semi minor axis of the ellipsoid in mm. 651 652 support : ART.ModuleSupport.Support 653 654 type : str 'Ellipsoidal Mirror'. 655 656 Methods 657 ------- 658 MirrorEllipsoidal.get_normal(Point) 659 660 MirrorEllipsoidal.get_centre() 661 662 MirrorEllipsoidal.get_grid3D(NbPoints) 663 664 """ 665 666 def __init__( 667 self, 668 Support, 669 SemiMajorAxis=None, 670 SemiMinorAxis=None, 671 OffAxisAngle=None, 672 f_object=None, 673 f_image=None, 674 IncidenceAngle=None, 675 Magnification=None, 676 DistanceObjectImage=None, 677 Eccentricity=None, 678 **kwargs, 679 ): 680 """ 681 Generate an ellipsoidal mirror with given parameters. 682 683 The angles are given in degrees but converted to radians for internal calculations. 684 You can specify the semi-major and semi-minor axes, the off-axis angle, the object and image focal distances or some other parameters. 685 The constructor uses a magic calculator to determine the missing parameters. 686 687 Parameters 688 ---------- 689 Support : TYPE 690 ART.ModuleSupport.Support. 691 SemiMajorAxis : float (optional) 692 Semi major axis of the ellipsoid in mm.. 693 SemiMinorAxis : float (optional) 694 Semi minor axis of the ellipsoid in mm.. 695 OffAxisAngle : float (optional) 696 Off-axis angle of the mirror in mm. Defined as the angle at the centre of the mirror between the two foci.. 697 f_object : float (optional) 698 Object focal distance in mm. 699 f_image : float (optional) 700 Image focal distance in mm. 701 IncidenceAngle : float (optional) 702 Angle of incidence in degrees. 703 Magnification : float (optional) 704 Magnification. 705 DistanceObjectImage : float (optional) 706 Distance between object and image in mm. 707 Eccentricity : float (optional) 708 Eccentricity of the ellipsoid. 709 """ 710 super().__init__() 711 self.type = "Ellipsoidal Mirror" 712 self.support = Support 713 714 if "Surface" in kwargs: 715 self.Surface = kwargs["Surface"] 716 717 parameter_calculator = EllipsoidalMirrorCalculator() 718 values, steps = parameter_calculator.calculate_values(verify_consistency=True, 719 f1=f_object, 720 f2=f_image, 721 a=SemiMajorAxis, 722 b=SemiMinorAxis, 723 offset_angle=OffAxisAngle, 724 incidence_angle=IncidenceAngle, 725 m=Magnification, 726 l=DistanceObjectImage, 727 e=Eccentricity) 728 self.a = values["a"] 729 self.b = values["b"] 730 self._offaxisangle = np.deg2rad(values["offset_angle"]) 731 self._f_object = values["f1"] 732 self._f_image = values["f2"] 733 734 self.centre_ref = self._get_centre_ref() 735 736 self.r0 = self.centre_ref 737 738 self.support_normal_ref = Vector([0, 0, 1.0]) 739 self.majoraxis_ref = Vector([1.0, 0, 0]) 740 741 self.f1_ref = Point([-self.a, 0,0]) 742 self.f2_ref = Point([ self.a, 0,0]) 743 self.towards_image_ref = (self.f2_ref - self.centre_ref).normalized() 744 self.towards_object_ref = (self.f1_ref - self.centre_ref).normalized() 745 self.centre_normal_ref = self.get_local_normal(self.centre_ref) 746 747 self.add_global_points("f1", "f2", "centre") 748 self.add_global_vectors("towards_image", "towards_object", "centre_normal", "support_normal", "majoraxis") 749 750 def _get_intersection(self, Ray): 751 """Return the intersection point between Ray and the ellipsoidal mirror surface.""" 752 ux, uy, uz = Ray.vector 753 xA, yA, zA = Ray.point 754 755 da = (uy**2 + uz**2) / self.b**2 + (ux / self.a) ** 2 756 db = 2 * ((uy * yA + uz * zA) / self.b**2 + (ux * xA) / self.a**2) 757 dc = (yA**2 + zA**2) / self.b**2 + (xA / self.a) ** 2 - 1 758 759 Solution = mgeo.SolverQuadratic(da, db, dc) 760 Solution = mgeo.KeepPositiveSolution(Solution) 761 762 ListPointIntersection = [] 763 C = self.get_centre_ref() 764 for t in Solution: 765 Intersect = Ray.vector * t + Ray.point 766 if Intersect[2] < 0 and Intersect - C in self.support: 767 ListPointIntersection.append(Intersect) 768 769 return _IntersectionRayMirror(Ray.point, ListPointIntersection) 770 771 def get_local_normal(self, Point): 772 """Return the normal unit vector on the ellipsoidal surface at point Point.""" 773 Gradient = Vector(np.zeros(3)) 774 775 Gradient[0] = -Point[0] / self.a**2 776 Gradient[1] = -Point[1] / self.b**2 777 Gradient[2] = -Point[2] / self.b**2 778 779 return Gradient.normalized() 780 781 def _get_centre_ref(self): 782 """Return 3D coordinates of the point on the mirror surface at the center of its support.""" 783 foci = 2 * np.sqrt(self.a**2 - self.b**2) # distance between the foci 784 h = -foci / 2 / np.tan(self._offaxisangle) 785 R = np.sqrt(foci**2 / 4 + h**2) 786 sign = 1 787 if math.isclose(self._offaxisangle, np.pi / 2): 788 h = 0 789 elif self._offaxisangle > np.pi / 2: 790 h = -h 791 sign = -1 792 a = 1 - self.a**2 / self.b**2 793 b = -2 * h 794 c = self.a**2 + h**2 - R**2 795 z = (-b + sign * np.sqrt(b**2 - 4 * a * c)) / (2 * a) 796 if math.isclose(z**2, self.b**2): 797 return np.array([0, 0, -self.b]) 798 x = self.a * np.sqrt(1 - z**2 / self.b**2) 799 centre = Point([x, 0, sign * z]) 800 return centre 801 802 def _zfunc(self, PointArray): 803 x = PointArray[:,0] 804 y = PointArray[:,1] 805 x-= self.r0[0] 806 y-= self.r0[1] 807 z = -np.sqrt(1 - (x / self.a)**2 - (y / self.b)**2) 808 z += self.r0[2] 809 return z 810 811 812# %% Cylindrical mirror definitions 813class MirrorCylindrical(Mirror): 814 """ 815 Cylindrical mirror surface with eqn. $y^2 + z^2 = R^2$, where $R$ is the radius. 816 817 Attributes 818 ---------- 819 radius : float 820 Radius of curvature. A postive value for concave mirror, a negative value for a convex mirror. 821 822 support : ART.ModuleSupport.Support 823 824 type : str 'Cylindrical Mirror'. 825 826 Methods 827 ------- 828 MirrorCylindrical.get_normal(Point) 829 830 MirrorCylindrical.get_centre() 831 832 MirrorCylindrical.get_grid3D(NbPoints) 833 """ 834 835 def __init__(self, Support, Radius, **kwargs): 836 """ 837 Construct a cylindrical mirror. 838 839 Parameters 840 ---------- 841 Radius : float 842 The radius of curvature in mm. A postive value for concave mirror, a negative value for a convex mirror. 843 844 Support : ART.ModuleSupport.Support 845 846 """ 847 super().__init__() 848 if Radius < 0: 849 self.type = "CylindricalCX Mirror" 850 self.curvature = Curvature.CONVEX 851 self.radius = -Radius 852 else: 853 self.type = "CylindricalCC Mirror" 854 self.curvature = Curvature.CONCAVE 855 self.radius = Radius 856 857 self.support = Support 858 859 if "Surface" in kwargs: 860 self.Surface = kwargs["Surface"] 861 862 self.r0 = Point([0.0, 0.0, -self.radius]) 863 864 self.support_normal_ref = Vector([0, 0, 1.0]) 865 self.majoraxis_ref = Vector([1.0, 0, 0]) 866 self.centre_ref = Point([0.0, 0.0, -self.radius]) 867 868 self.add_global_points("centre") 869 self.add_global_vectors("support_normal", "majoraxis") 870 871 def _get_intersection(self, Ray): 872 """Return the intersection point between the Ray and the cylinder.""" 873 uy = Ray.vector[1] 874 uz = Ray.vector[2] 875 yA = Ray.point[1] 876 zA = Ray.point[2] 877 878 a = uy**2 + uz**2 879 b = 2 * (uy * yA + uz * zA) 880 c = yA**2 + zA**2 - self.radius**2 881 882 Solution = mgeo.SolverQuadratic(a, b, c) 883 Solution = mgeo.KeepPositiveSolution(Solution) 884 885 ListPointIntersection = [] 886 for t in Solution: 887 Intersect = Ray.vector * t + Ray.point 888 if Intersect[2] < 0 and Intersect in self.support: 889 ListPointIntersection.append(Intersect) 890 891 return _IntersectionRayMirror(Ray.point, ListPointIntersection) 892 893 def get_local_normal(self, Point): 894 """Return the normal unit vector on the cylinder surface at point P.""" 895 Gradient = Vector([0, -Point[1], -Point[2]]) 896 return Gradient.normalized() 897 898 def get_centre(self): 899 """Return 3D coordinates of the point on the mirror surface at the center of its support.""" 900 return Point([0, 0, -self.radius]) 901 902 def _zfunc(self, PointArray): 903 y = PointArray[:,1] 904 return -np.sqrt(self.radius**2 - y**2) 905 906 907# %% Grazing parabolic mirror definitions 908# A grazing parabola is no different from a parabola, it's just a question of what we consider to be the support. 909# 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. 910 911class GrazingParabola(Mirror): 912 r""" 913 A parabolic mirror with a support parallel to the surface at the optical center. 914 Needs to be completed, currently not functional. 915 TODO 916 """ 917 918 def __init__(self, Support, 919 FocalEffective: float=None, 920 OffAxisAngle: float = None, 921 FocalParent: float = None, 922 RadiusParent: float = None, 923 OffAxisDistance: float = None, 924 MoreThan90: bool = None, 925 **kwargs): 926 """ 927 Initialise a Parabolic mirror. 928 929 Parameters 930 ---------- 931 FocalEffective : float 932 Effective focal length of the parabola in mm. 933 934 OffAxisAngle : float 935 Off-axis angle *in degrees* of the parabola. 936 937 Support : ART.ModuleSupport.Support 938 939 """ 940 super().__init__() 941 self.curvature = Curvature.CONCAVE 942 self.support = Support 943 self.type = "Grazing Parabolic Mirror" 944 945 if "Surface" in kwargs: 946 self.Surface = kwargs["Surface"] 947 948 parameter_calculator = OAP_calculator() 949 values,steps = parameter_calculator.calculate_values(verify_consistency=True, 950 fs=FocalEffective, 951 theta=OffAxisAngle, 952 fp=FocalParent, 953 Rc=RadiusParent, 954 OAD=OffAxisDistance, 955 more_than_90=MoreThan90) 956 self._offaxisangle = values["theta"] 957 self._feff = values["fs"] 958 self._p = values["p"] 959 self._offaxisdistance = values["OAD"] 960 self._fparent = values["fp"] 961 self._rparent = values["Rc"] 962 963 #self.r0 = Point([self._offaxisdistance, 0, self._offaxisdistance**2 / 2 / self._p]) 964 self.r0 = Point([0, 0, 0]) 965 966 self.centre_ref = Point([0, 0, 0]) 967 self.support_normal_ref = Vector([0, 0, 1.0]) 968 self.majoraxis_ref = Vector([1.0, 0, 0]) 969 970 focus_ref = Point([np.sqrt(self._feff**2 - self._offaxisdistance**2), 0, self._offaxisdistance]) # frame where x is as collimated direction 971 # We need to rotate it in the trigonometric direction by the 90°-theta/2 972 angle = np.pi/2 - self._offaxisangle/2 973 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), ]) 974 975 self.towards_focusing_ref = (self.focus_ref-self.centre_ref).normalized() 976 towards_collimated_ref = Vector([-1.0, 0, 0]) # Again, we need to rotate it by the same angle 977 self.collimated_ref = Point([-np.cos(angle), 0, np.sin(angle)]) 978 979 self.add_global_points("focus", "centre") 980 self.add_global_vectors("towards_focusing", "towards_collimated", "support_normal", "majoraxis") 981 982 def _get_intersection(self, Ray): 983 """ 984 Return the intersection point between the ray and the parabola. 985 """ 986 pass 987 988 def get_local_normal(self, Point): 989 """ 990 Return the normal unit vector on the paraboloid surface at point Point. 991 """ 992 pass 993 994 def _zfunc(self, PointArray): 995 pass 996 997 998# %% Reflections on mirrors, is it useful to keep this? 999def _ReflectionMirrorRay(Mirror, PointMirror, Ray): 1000 """ 1001 Return the reflected ray according to the law of reflection. 1002 1003 Parameters 1004 ---------- 1005 Mirror : Mirror-objectS 1006 1007 PointMirror : np.ndarray 1008 Point of reflection on the mirror surface. 1009 1010 Ray : Ray-object 1011 1012 """ 1013 PointRay = Ray.point 1014 VectorRay = Ray.vector 1015 NormalMirror = Mirror.get_local_normal(PointMirror) 1016 1017 #VectorRayReflected = mgeo.SymmetricalVector(-VectorRay, NormalMirror) 1018 VectorRayReflected = VectorRay- 2*NormalMirror*np.dot(VectorRay,NormalMirror) # Is it any better than SymmetricalVector? 1019 1020 RayReflected = Ray.copy_ray() 1021 RayReflected.point = PointMirror 1022 RayReflected.vector = VectorRayReflected 1023 RayReflected.incidence = mgeo.AngleBetweenTwoVectors( 1024 -VectorRay, NormalMirror 1025 ) 1026 RayReflected.path = Ray.path + (np.linalg.norm(PointMirror - PointRay),) 1027 1028 return RayReflected 1029 1030 1031def ReflectionMirrorRayList(Mirror, ListRay, IgnoreDefects=False): 1032 """ 1033 Return the the reflected rays according to the law of reflection for the list of incident rays ListRay. 1034 1035 Rays that do not hit the support are not further propagated. 1036 1037 Updates the reflected rays' incidence angle and path. 1038 1039 Parameters 1040 ---------- 1041 Mirror : Mirror-object 1042 1043 ListRay : list[Ray-object] 1044 1045 """ 1046 Deformed = type(Mirror) == DeformedMirror 1047 ListRayReflected = [] 1048 for k in ListRay: 1049 PointMirror = Mirror._get_intersection(k) 1050 1051 if PointMirror is not None: 1052 if Deformed and IgnoreDefects: 1053 M = Mirror.Mirror 1054 else: 1055 M = Mirror 1056 RayReflected = _ReflectionMirrorRay(M, PointMirror, k) 1057 ListRayReflected.append(RayReflected) 1058 return ListRayReflected 1059 1060 1061# %% Deformed mirror definitions, need to be reworked, maybe implemented as coatings? 1062class DeformedMirror(Mirror): 1063 def __init__(self, Mirror, DeformationList): 1064 self.Mirror = Mirror 1065 self.DeformationList = DeformationList 1066 self.type = Mirror.type 1067 self.support = self.Mirror.support 1068 1069 def get_local_normal(self, PointMirror): 1070 base_normal = self.Mirror.get_normal(PointMirror) 1071 C = self.get_centre() 1072 defects_normals = [ 1073 d.get_normal(PointMirror - C) for d in self.DeformationList 1074 ] 1075 for i in defects_normals: 1076 base_normal = mgeo.normal_add(base_normal, i) 1077 base_normal /= np.linalg.norm(base_normal) 1078 return base_normal 1079 1080 def get_centre(self): 1081 return self.Mirror.get_centre() 1082 1083 def get_grid3D(self, NbPoint, **kwargs): 1084 return self.Mirror.get_grid3D(NbPoint, **kwargs) 1085 1086 def _get_intersection(self, Ray): 1087 Intersect = self.Mirror._get_intersection(Ray) 1088 if Intersect is not None: 1089 h = sum( 1090 D.get_offset(Intersect - self.get_centre()) 1091 for D in self.DeformationList 1092 ) 1093 alpha = mgeo.AngleBetweenTwoVectors( 1094 -Ray.vector, self.Mirror.get_normal(Intersect) 1095 ) 1096 Intersect -= Ray.vector * h / np.cos(alpha) 1097 return Intersect
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.
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.
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.
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.
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.
654def DrawAsphericity(Mirror, Npoints=1000): 655 """ 656 This function displays a map of the asphericity of the mirror. 657 It's a scatter plot of the points of the mirror surface, with the color representing the distance to the closest sphere. 658 The closest sphere is calculated by the function get_closest_sphere, so least square method. 659 660 Parameters 661 ---------- 662 Mirror : Mirror 663 The mirror to analyse. 664 665 Npoints : int, optional 666 The number of points to sample on the mirror surface. The default is 1000. 667 668 Returns 669 ------- 670 fig : Figure 671 The figure of the plot. 672 """ 673 plt.ion() 674 fig = plt.figure() 675 ax = Mirror.support._ContourSupport(fig) 676 center, radius = man.GetClosestSphere(Mirror, Npoints) 677 Points = mpm.sample_support(Mirror.support, Npoints=1000) 678 Points += Mirror.r0[:2] 679 Z = Mirror._zfunc(Points) 680 Points = mgeo.PointArray([Points[:, 0], Points[:, 1], Z]).T 681 X, Y = Points[:, 0] - Mirror.r0[0], Points[:, 1] - Mirror.r0[1] 682 Points_centered = Points - center 683 Distance = np.linalg.norm(Points_centered, axis=1) - radius 684 Distance*=1e3 # To convert to µm 685 p = plt.scatter(X, Y, c=Distance, s=15) 686 divider = man.make_axes_locatable(ax) 687 cax = divider.append_axes("right", size="5%", pad=0.05) 688 cbar = fig.colorbar(p, cax=cax) 689 cbar.set_label("Distance to closest sphere (µm)") 690 ax.set_xlabel("x (mm)") 691 ax.set_ylabel("y (mm)") 692 plt.title("Asphericity map", loc="right") 693 plt.tight_layout() 694 695 bbox = ax.get_position() 696 bbox.set_points(bbox.get_points() - np.array([[0.01, 0], [0.01, 0]])) 697 ax.set_position(bbox) 698 plt.show() 699 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)
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
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.
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 for t in Solution: 289 Intersect = Ray.vector * t + Ray.point 290 if Intersect - self.r0 in self.support: 291 return Intersect, True 292 293 return Ray.point, False 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)
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)
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
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()
Return the normal unit vector on the spherical surface at point Point.
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
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
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)
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")
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
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()
Return the normal unit vector on the paraboloid surface at point Point.
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()
Vectorised version of the normal calculation.
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
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 )
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)
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")
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
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()
Return the normal unit vector on the toroidal surface at point Point in the reference frame of the mirror
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
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
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.
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
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)
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")
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.
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()
Return the normal unit vector on the ellipsoidal surface at point Point.
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
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)
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)
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")
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
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()
Return the normal unit vector on the cylinder surface at point P.
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])
Return 3D coordinates of the point on the mirror surface at the center of its support.
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
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
A parabolic mirror with a support parallel to the surface at the optical center. Needs to be completed, currently not functional. TODO
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")
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
989 def get_local_normal(self, Point): 990 """ 991 Return the normal unit vector on the paraboloid surface at point Point. 992 """ 993 pass
Return the normal unit vector on the paraboloid surface at point Point.
Inherited Members
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
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]
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
Abstract base class for mirrors.
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
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.