ModuleOpticalChain
Provides the class OpticalChain, which represents the whole optical setup to be simulated, from the bundle of source-Rays through the successive OpticalElements.
Created in Jan 2022
@author: Stefan Haessler
1""" 2Provides the class OpticalChain, which represents the whole optical setup to be simulated, 3from the bundle of source-[Rays](ModuleOpticalRay.html) through the successive [OpticalElements](ModuleOpticalElement.html). 4 5 6 7 8 9Created in Jan 2022 10 11@author: Stefan Haessler 12""" 13# %% Modules 14import copy 15import numpy as np 16import logging 17 18import ARTcore.ModuleProcessing as mp 19import ARTcore.ModuleGeometry as mgeo 20import ARTcore.ModuleOpticalRay as mray 21import ARTcore.ModuleOpticalElement as moe 22import ARTcore.ModuleDetector as mdet 23import ARTcore.ModuleSource as msource 24 25logger = logging.getLogger(__name__) 26 27# %% 28class OpticalChain: 29 """ 30 The OpticalChain represents the whole optical setup to be simulated: 31 Its main attributes are a list source-[Rays](ModuleOpticalRay.html) and 32 a list of successive [OpticalElements](ModuleOpticalElement.html). 33 34 The method OpticalChain.get_output_rays() returns an associated list of lists of 35 [Rays](ModuleOpticalRay.html), each calculated by ray-tracing from one 36 [OpticalElement](ModuleOpticalElement.html) to the next. 37 So OpticalChain.get_output_rays()[i] is the bundle of [Rays](ModuleOpticalRay.html) *after* 38 optical_elements[i]. 39 40 The string "description" can contain a short description of the optical setup, or similar notes. 41 42 The OpticalChain can be visualized quickly with the method OpticalChain.quickshow(), 43 and more nicely with OpticalChain.render(). 44 45 The class also provides methods for (mis-)alignment of the source-[Ray](ModuleOpticalRay.html)-bundle and the 46 [OpticalElements](ModuleOpticalElement.html). 47 48 Attributes 49 ---------- 50 source_rays : list[mray.Ray] 51 List of source rays, which are to be traced. 52 53 optical_elements : list[moe.OpticalElement] 54 List of successive optical elements. 55 56 detector: mdet.Detector (optional) 57 The detector (or list of detectors) to analyse the results. 58 59 description : str 60 A string to describe the optical setup. 61 Methods 62 ---------- 63 copy_chain() 64 65 get_output_rays() 66 67 render() 68 69 ---------- 70 71 shift_source(axis, distance) 72 73 tilt_source(self, axis, angle) 74 75 ---------- 76 77 rotate_OE(OEindx, axis, angle) 78 79 shift_OE(OEindx, axis, distance) 80 81 """ 82 83 def __init__(self, source_rays, optical_elements, detectors, description=""): 84 """ 85 Parameters 86 ---------- 87 source_rays : list[mray.Ray] 88 List of source rays, which are to be traced. 89 90 optical_elements : list[moe.OpticalElement] 91 List of successive optical elements. 92 93 detector: mdet.Detector (optional) 94 The detector (or list of detectors) to analyse the results. 95 96 description : str, optional 97 A string to describe the optical setup. Defaults to ''. 98 99 """ 100 self.source_rays = copy.deepcopy(source_rays) 101 # deepcopy so this object doesn't get changed when the global source_rays change "outside" 102 self.optical_elements = copy.deepcopy(optical_elements) 103 # deepcopy so this object doesn't get changed when the global optical_elements changes "outside" 104 self.detectors = detectors 105 if isinstance(detectors, mdet.Detector): 106 self.detectors = {"Focus": detectors} 107 self.description = description 108 self._output_rays = None 109 self._last_source_rays_hash = None 110 self._last_optical_elements_hash = None 111 112 def __repr__(self): 113 pretty_str = "Optical setup [OpticalChain]:\n" 114 pretty_str += f" - Description: {self.description}\n" if self.description else "Description: Not provided.\n" 115 pretty_str += " - Contains the following elements:\n" 116 pretty_str += f" - Source with {len(self.source_rays)} rays at coordinate origin\n" 117 prev_pos = np.zeros(3) 118 for i, element in enumerate(self.optical_elements): 119 dist = (element.position - prev_pos).norm 120 if i == 0: 121 prev = "source" 122 else: 123 prev = f"element {i-1}" 124 pretty_str += f" - Element {i}: {element.type} at distance {round(dist)} from {prev}\n" 125 prev_pos = element.position 126 for i in self.detectors.keys(): 127 detector = self.detectors[i] 128 pretty_str += f' - Detector "{i}": {detector.type} at distance {round(detector.distance,3)} from element {detector.index % len(self.optical_elements)}\n' 129 return pretty_str 130 131 def __getitem__(self, i): 132 return self.optical_elements[i] 133 134 def __len__(self): 135 return len(self.optical_elements) 136 137 @property 138 def source_rays(self): 139 return self._source_rays 140 141 @source_rays.setter 142 def source_rays(self, source_rays): 143 if type(source_rays) == mray.RayList: 144 self._source_rays = source_rays 145 else: 146 raise TypeError("Source_rays must be a RayList object.") 147 148 @property 149 def optical_elements(self): 150 return self._optical_elements 151 152 @optical_elements.setter 153 def optical_elements(self, optical_elements): 154 if type(optical_elements) == list and all(isinstance(x, moe.OpticalElement) for x in optical_elements): 155 self._optical_elements = optical_elements 156 else: 157 raise TypeError("Optical_elements must be list of OpticalElement-objects.") 158 159 # %% METHODS ################################################################## 160 161 def __copy__(self): 162 """Return another optical chain with the same source, optical elements and description-string as this one.""" 163 return OpticalChain(self.source_rays, self.optical_elements, self.detectors, self.description) 164 165 def __deepcopy__(self, memo): 166 """Return another optical chain with the same source, optical elements and description-string as this one.""" 167 return OpticalChain(copy.deepcopy(self.source_rays), copy.deepcopy(self.optical_elements), memo+"\n"+copy.copy(self.description)) 168 169 def get_input_rays(self): 170 """ 171 Returns the list of source rays. 172 """ 173 return self.source_rays 174 175 def get_output_rays(self, force=False, **kwargs): 176 """ 177 Returns the list of (lists of) output rays, calculate them if this hasn't been done yet, 178 or if the source-ray-bundle or anything about the optical elements has changed. 179 180 This is the user-facing method to perform the ray-tracing calculation. 181 """ 182 current_source_rays_hash = hash(self.source_rays) 183 current_optical_elements_hash = mp._hash_list_of_objects(self.optical_elements) 184 if (current_source_rays_hash != self._last_source_rays_hash) or (current_optical_elements_hash != self._last_optical_elements_hash) or force: 185 print("...ray-tracing...", end="", flush=True) 186 self._output_rays = mp.RayTracingCalculation(self.source_rays, self.optical_elements, **kwargs) 187 print( 188 "\r\033[K", end="", flush=True 189 ) # move to beginning of the line with \r and then delete the whole line with \033[K 190 self._last_source_rays_hash = current_source_rays_hash 191 self._last_optical_elements_hash = current_optical_elements_hash 192 193 return self._output_rays 194 195 def get2dPoints(self, Detector = "Focus"): 196 """ 197 Returns the 2D points of the detected rays on the detector. 198 """ 199 if isinstance(Detector, str): 200 if Detector in self.detectors.keys(): 201 Detector = self.detectors[Detector] 202 else: 203 raise ValueError('The "Detector"-argument must be a name of an existing detector or a Detector-object.') 204 return Detector.get_2D_points(self.get_output_rays()[Detector.index])[0] 205 206 def get3dPoints(self, Detector = "Focus"): 207 """ 208 Returns the 3D points of the detected rays on the detector. 209 """ 210 if isinstance(Detector, str): 211 if Detector in self.detectors.keys(): 212 Detector = self.detectors[Detector] 213 else: 214 raise ValueError('The "Detector"-argument must be a name of an existing detector or a Detector-object.') 215 return Detector.get_3D_points(self.get_output_rays()[Detector.index]) 216 217 def getDelays(self, Detector = "Focus"): 218 """ 219 Returns the delays of the detected rays on the detector. 220 """ 221 if isinstance(Detector, str): 222 if Detector in self.detectors.keys(): 223 Detector = self.detectors[Detector] 224 else: 225 raise ValueError('The "Detector"-argument must be a name of an existing detector or a Detector-object.') 226 return Detector.get_Delays(self.get_output_rays()[Detector.index]) 227 228 # %% methods to (mis-)align the optical chain; just uses the corresponding methods of the OpticalElement class... 229 230 def shift_source(self, axis: str|np.ndarray, distance: float): 231 """ 232 Shift source ray bundle by distance (in mm) along the 'axis' specified as 233 a lab-frame vector (numpy-array of length 3) or as one of the strings 234 "vert", "horiz", or "random". 235 236 In the latter case, the reference is the incidence plane of the first 237 non-normal-incidence mirror after the source. If there is none, you will 238 be asked to rather specify the axis as a 3D-numpy-array. 239 240 axis = "vert" means the source position is shifted along the axis perpendicular 241 to that incidence plane, i.e. "vertically" away from the former incidence plane.. 242 243 axis = "horiz" means the source direciton is translated along thr axis in that 244 incidence plane and perpendicular to the current source direction, 245 i.e. "horizontally" in the incidence plane, but retaining the same distance 246 of source and first optical element. 247 248 axis = "random" means the the source direction shifted in a random direction 249 within in the plane perpendicular to the current source direction, 250 e.g. simulating a fluctuation of hte transverse source position. 251 252 Parameters 253 ---------- 254 axis : np.ndarray or str 255 Shift axis, specified either as a 3D lab-frame vector or as one 256 of the strings "vert", "horiz", or "random". 257 258 distance : float 259 Shift distance in mm. 260 261 Returns 262 ------- 263 Nothing, just modifies the property 'source_rays'. 264 """ 265 if type(distance) not in [int, float, np.float64]: 266 raise ValueError('The "distance"-argument must be an int or float number.') 267 268 central_ray_vector = mp.FindCentralRay(self.source_rays).vector 269 mirror_indcs = [i for i, OE in enumerate(self.optical_elements) if "Mirror" in OE.type.type] 270 271 OEnormal = None 272 for i in mirror_indcs: 273 ith_OEnormal = self.optical_elements[i].normal 274 if np.linalg.norm(np.cross(central_ray_vector, ith_OEnormal)) > 1e-10: 275 OEnormal = ith_OEnormal 276 break 277 if OEnormal is None: 278 raise Exception( 279 "There doesn't seem to be a non-normal-incidence mirror in this optical chain, \ 280 so you should rather give 'axis' as a numpy-array of length 3." 281 ) 282 283 if type(axis) == np.ndarray and len(axis) == 3: 284 translation_vector = axis 285 else: 286 perp_axis = np.cross(central_ray_vector, OEnormal) 287 horiz_axis = np.cross(perp_axis, central_ray_vector) 288 289 if axis == "vert": 290 translation_vector = perp_axis 291 elif axis == "horiz": 292 translation_vector = horiz_axis 293 elif axis == "random": 294 translation_vector = ( 295 np.random.uniform(low=-1, high=1, size=1) * perp_axis 296 + np.random.uniform(low=-1, high=1, size=1) * horiz_axis 297 ) 298 else: 299 raise ValueError( 300 'The shift direction must be specified by "axis" as one of ["vert", "horiz", "random"].' 301 ) 302 303 self.source_rays = mgeo.TranslationRayList(self.source_rays, distance * mgeo.Normalize(translation_vector)) 304 305 def tilt_source(self, axis: str|np.ndarray, angle: float): 306 """ 307 Rotate source ray bundle by angle around an axis, specified as 308 a lab-frame vector (numpy-array of length 3) or as one of the strings 309 "in_plane", "out_plane" or "random" direction. 310 311 In the latter case, the function considers the incidence plane of the first 312 non-normal-incidence mirror after the source. If there is none, you will 313 be asked to rather specify the axis as a 3D-numpy-array. 314 315 axis = "in_plane" means the source direction is rotated about an axis 316 perpendicular to that incidence plane, which tilts the source 317 "horizontally" in the same plane. 318 319 axis = "out_plane" means the source direciton is rotated about an axis 320 in that incidence plane and perpendicular to the current source direction, 321 which tilts the source "vertically" out of the former incidence plane. 322 323 axis = "random" means the the source direction is tilted in a random direction, 324 e.g. simulating a beam pointing fluctuation. 325 326 Attention, "angle" is given in deg, so as to remain consitent with the 327 conventions of other functions, although pointing is mostly talked about 328 in mrad instead. 329 330 Parameters 331 ---------- 332 axis : np.ndarray or str 333 Shift axis, specified either as a 3D lab-frame vector or as one 334 of the strings "in_plane", "out_plane", or "random". 335 336 angle : float 337 Rotation angle in degree. 338 339 Returns 340 ------- 341 Nothing, just modifies the property 'source_rays'. 342 """ 343 if type(angle) not in [int, float, np.float64]: 344 raise ValueError('The "angle"-argument must be an int or float number.') 345 346 central_ray_vector = mp.FindCentralRay(self.source_rays).vector 347 mirror_indcs = [i for i, OE in enumerate(self.optical_elements) if "Mirror" in OE.type.type] 348 349 OEnormal = None 350 for i in mirror_indcs: 351 ith_OEnormal = self.optical_elements[i].normal 352 if np.linalg.norm(np.cross(central_ray_vector, ith_OEnormal)) > 1e-10: 353 OEnormal = ith_OEnormal 354 break 355 if OEnormal is None: 356 raise Exception( 357 "There doesn't seem to be a non-normal-incidence mirror in this optical chain, \ 358 so you should rather give 'axis' as a numpy-array of length 3." 359 ) 360 361 if type(axis) == np.ndarray and len(axis) == 3: 362 rot_axis = axis 363 else: 364 rot_axis_in = np.cross(central_ray_vector, OEnormal) 365 rot_axis_out = np.cross(rot_axis_in, central_ray_vector) 366 if axis == "in_plane": 367 rot_axis = rot_axis_in 368 elif axis == "out_plane": 369 rot_axis = rot_axis_out 370 elif axis == "random": 371 rot_axis = ( 372 np.random.uniform(low=-1, high=1, size=1) * rot_axis_in 373 + np.random.uniform(low=-1, high=1, size=1) * rot_axis_out 374 ) 375 else: 376 raise ValueError( 377 'The tilt axis must be specified by as one of ["in_plane", "out_plane", "random"] or as a numpy-array of length 3.' 378 ) 379 380 self.source_rays = mgeo.RotationAroundAxisRayList(self.source_rays, rot_axis, np.deg2rad(angle)) 381 382 def partial_realign(self, OEstart: int, OEstop: int, DistanceList: list, IncidenceAngleList, IncidencePlaneAngleList): 383 """ 384 This as-of-yet not-implemented method will realign only the parts of the optical chain that are between the elements 385 OEstart and OEstop (both included). 386 """ 387 OpticsList = [i.type for i in self.optical_elements[OEstart:OEstop]] 388 outray = self.get_output_rays()[OEstart-1][0] 389 print(outray) 390 new_elements = mp.OEPlacement(OpticsList, DistanceList, IncidenceAngleList, IncidencePlaneAngleList, outray) 391 self.optical_elements[OEstart:OEstop] = new_elements 392 393 # %% 394 def rotate_OE(self, OEindx: int, ref: str, axis: str, angle: float): 395 """ 396 Rotate the optical element OpticalChain.optical_elements[OEindx] with an axis defined relative to the master ray. 397 The axis can either be relative to the incoming master ray or the outgoing one. This is defined by the "ref" variable 398 that takes either of the two values: 399 - "in" 400 - "out" 401 In either case the "axis" can take these values: 402 - "roll": Rotation about the master ray. Looking in the same direction as light propagation, positive is counterclockwise 403 - "pitch": Rotation about the vector in the incidence plane that is normal to the master ray. 404 - "yaw": Rotation about the vector normal to the incidence plane and to the master ray. 405 406 Parameters 407 ---------- 408 OEindx : int 409 Index of the optical element to modify out of OpticalChain.optical_elements. 410 411 ref : str 412 Reference ray used to define the axes of rotation. Can be either: 413 "in" or "out" or "local_normal" (in+out) 414 415 axis : str 416 Rotation axis, specified as one of the strings 417 "pitch", "roll", "yaw" 418 419 angle : float 420 Rotation angle in degree. 421 422 Returns 423 ------- 424 Nothing, just modifies OpticalChain.optical_elements[OEindx]. 425 """ 426 if abs(OEindx) > len(self.optical_elements): 427 raise ValueError( 428 'The "OEnumber"-argument is out of range compared to the length of OpticalChain.optical_elements.' 429 ) 430 if type(angle) not in [int, float, np.float64]: 431 raise ValueError('The "angle"-argument must be an int or float number.') 432 MasterRay = [mp.FindCentralRay(self.source_rays)] 433 TracedRay = mp.RayTracingCalculation(MasterRay, self.optical_elements) 434 TracedRay = [MasterRay] + TracedRay 435 match ref: 436 case "out": 437 RefVec = mgeo.Normalize(TracedRay[OEindx+1][0].vector) 438 case "in": 439 RefVec = mgeo.Normalize(TracedRay[OEindx][0].vector) 440 case "localnormal": 441 In = mgeo.Normalize(TracedRay[OEindx][0].vector) 442 Out = mgeo.Normalize(TracedRay[OEindx+1][0].vector) 443 RefVec = mgeo.Normalize(Out-In) 444 case _: 445 raise ValueError('The "ref"-argument must be a string out of ["in", "out", "localnormal].') 446 if 1 - np.dot(self[OEindx].normal,RefVec) > 1e-2: 447 match axis: 448 case "pitch": 449 self[OEindx].normal = mgeo.RotationAroundAxis(self[OEindx].majoraxis,np.deg2rad(angle),self[OEindx].normal) 450 case "roll": 451 self[OEindx].normal = mgeo.RotationAroundAxis(RefVec,np.deg2rad(angle),self[OEindx].normal) 452 case "yaw": 453 self[OEindx].normal = mgeo.RotationAroundAxis(mgeo.Normalize(np.cross(RefVec, self[OEindx].majoraxis)),np.deg2rad(angle),self[OEindx].normal) 454 case _: 455 raise ValueError('The "axis"-argument must be a string out of ["pitch", "roll", "yaw"].') 456 else: 457 #If the normal vector is aligned with the ray 458 match axis: 459 case "pitch": 460 self[OEindx].normal = mgeo.RotationAroundAxis(self[OEindx].majoraxis,np.deg2rad(angle),self[OEindx].normal) 461 case "roll": 462 self[OEindx].majoraxis = mgeo.RotationAroundAxis(self[OEindx].normal,np.deg2rad(angle),self[OEindx].majoraxis) 463 case "yaw": 464 self[OEindx].normal = mgeo.RotationAroundAxis(mgeo.Normalize(np.cross(RefVec, self[OEindx].majoraxis)),np.deg2rad(angle),self[OEindx].normal) 465 case _: 466 raise ValueError('The "axis"-argument must be a string out of ["pitch", "roll", "yaw"].') 467 468 def shift_OE(self, OEindx: int, ref: str, axis: str, distance: float): 469 """ 470 Shift the optical element OpticalChain.optical_elements[OEindx] along 471 axis referenced to the master ray. 472 473 The axis can either be relative to the incoming master ray or the outgoing one. This is defined by the "ref" variable 474 that takes either of the two values: 475 - "in" 476 - "out" 477 In either case the "axis" can take these values: 478 - "along": Translation along the master ray. 479 - "in_plane": Translation along the vector in the incidence plane that is normal to the master ray. 480 - "out_plane": Translation along the vector normal to the incidence plane and to the master ray. 481 482 Parameters 483 ---------- 484 OEindx : int 485 Index of the optical element to modify out of OpticalChain.optical_elements. 486 487 ref : str 488 Reference ray used to define the axes of rotation. Can be either: 489 "in" or "out". 490 491 axis : str 492 Translation axis, specified as one of the strings 493 "along", "in_plane", "out_plane". 494 495 distance : float 496 Rotation angle in degree. 497 498 Returns 499 ------- 500 Nothing, just modifies OpticalChain.optical_elements[OEindx]. 501 """ 502 if abs(OEindx) >= len(self.optical_elements): 503 raise ValueError( 504 'The "OEnumber"-argument is out of range compared to the length of OpticalChain.optical_elements.' 505 ) 506 if type(distance) not in [int, float, np.float64]: 507 raise ValueError('The "dist"-argument must be an int or float number.') 508 509 MasterRay = [mp.FindCentralRay(self.source_rays)] 510 TracedRay = mp.RayTracingCalculation(MasterRay, self.optical_elements) 511 TracedRay = [MasterRay] + TracedRay 512 match ref: 513 case "out": 514 RefVec = TracedRay[OEindx+1][0].vector 515 case "in": 516 RefVec = TracedRay[OEindx][0].vector 517 case _: 518 raise ValueError('The "ref"-argument must be a string out of ["in", "out"].') 519 match axis: 520 case "along": 521 self[OEindx].position = self[OEindx].position + distance * mgeo.Normalize(RefVec) 522 case "in_plane": 523 self[OEindx].position = self[OEindx].position + distance * mgeo.Normalize(np.cross(RefVec, np.cross(RefVec, self[OEindx].normal))) 524 case "out_plane": 525 self[OEindx].position = self[OEindx].position + distance * mgeo.Normalize(np.cross(RefVec, self[OEindx].normal)) 526 case _: 527 raise ValueError('The "axis"-argument must be a string out of ["along", "in_plane", "out_plane"].') 528 529 # some function that randomly misalings one, or several or all ?
logger =
<Logger ARTcore.ModuleOpticalChain (WARNING)>
class
OpticalChain:
29class OpticalChain: 30 """ 31 The OpticalChain represents the whole optical setup to be simulated: 32 Its main attributes are a list source-[Rays](ModuleOpticalRay.html) and 33 a list of successive [OpticalElements](ModuleOpticalElement.html). 34 35 The method OpticalChain.get_output_rays() returns an associated list of lists of 36 [Rays](ModuleOpticalRay.html), each calculated by ray-tracing from one 37 [OpticalElement](ModuleOpticalElement.html) to the next. 38 So OpticalChain.get_output_rays()[i] is the bundle of [Rays](ModuleOpticalRay.html) *after* 39 optical_elements[i]. 40 41 The string "description" can contain a short description of the optical setup, or similar notes. 42 43 The OpticalChain can be visualized quickly with the method OpticalChain.quickshow(), 44 and more nicely with OpticalChain.render(). 45 46 The class also provides methods for (mis-)alignment of the source-[Ray](ModuleOpticalRay.html)-bundle and the 47 [OpticalElements](ModuleOpticalElement.html). 48 49 Attributes 50 ---------- 51 source_rays : list[mray.Ray] 52 List of source rays, which are to be traced. 53 54 optical_elements : list[moe.OpticalElement] 55 List of successive optical elements. 56 57 detector: mdet.Detector (optional) 58 The detector (or list of detectors) to analyse the results. 59 60 description : str 61 A string to describe the optical setup. 62 Methods 63 ---------- 64 copy_chain() 65 66 get_output_rays() 67 68 render() 69 70 ---------- 71 72 shift_source(axis, distance) 73 74 tilt_source(self, axis, angle) 75 76 ---------- 77 78 rotate_OE(OEindx, axis, angle) 79 80 shift_OE(OEindx, axis, distance) 81 82 """ 83 84 def __init__(self, source_rays, optical_elements, detectors, description=""): 85 """ 86 Parameters 87 ---------- 88 source_rays : list[mray.Ray] 89 List of source rays, which are to be traced. 90 91 optical_elements : list[moe.OpticalElement] 92 List of successive optical elements. 93 94 detector: mdet.Detector (optional) 95 The detector (or list of detectors) to analyse the results. 96 97 description : str, optional 98 A string to describe the optical setup. Defaults to ''. 99 100 """ 101 self.source_rays = copy.deepcopy(source_rays) 102 # deepcopy so this object doesn't get changed when the global source_rays change "outside" 103 self.optical_elements = copy.deepcopy(optical_elements) 104 # deepcopy so this object doesn't get changed when the global optical_elements changes "outside" 105 self.detectors = detectors 106 if isinstance(detectors, mdet.Detector): 107 self.detectors = {"Focus": detectors} 108 self.description = description 109 self._output_rays = None 110 self._last_source_rays_hash = None 111 self._last_optical_elements_hash = None 112 113 def __repr__(self): 114 pretty_str = "Optical setup [OpticalChain]:\n" 115 pretty_str += f" - Description: {self.description}\n" if self.description else "Description: Not provided.\n" 116 pretty_str += " - Contains the following elements:\n" 117 pretty_str += f" - Source with {len(self.source_rays)} rays at coordinate origin\n" 118 prev_pos = np.zeros(3) 119 for i, element in enumerate(self.optical_elements): 120 dist = (element.position - prev_pos).norm 121 if i == 0: 122 prev = "source" 123 else: 124 prev = f"element {i-1}" 125 pretty_str += f" - Element {i}: {element.type} at distance {round(dist)} from {prev}\n" 126 prev_pos = element.position 127 for i in self.detectors.keys(): 128 detector = self.detectors[i] 129 pretty_str += f' - Detector "{i}": {detector.type} at distance {round(detector.distance,3)} from element {detector.index % len(self.optical_elements)}\n' 130 return pretty_str 131 132 def __getitem__(self, i): 133 return self.optical_elements[i] 134 135 def __len__(self): 136 return len(self.optical_elements) 137 138 @property 139 def source_rays(self): 140 return self._source_rays 141 142 @source_rays.setter 143 def source_rays(self, source_rays): 144 if type(source_rays) == mray.RayList: 145 self._source_rays = source_rays 146 else: 147 raise TypeError("Source_rays must be a RayList object.") 148 149 @property 150 def optical_elements(self): 151 return self._optical_elements 152 153 @optical_elements.setter 154 def optical_elements(self, optical_elements): 155 if type(optical_elements) == list and all(isinstance(x, moe.OpticalElement) for x in optical_elements): 156 self._optical_elements = optical_elements 157 else: 158 raise TypeError("Optical_elements must be list of OpticalElement-objects.") 159 160 # %% METHODS ################################################################## 161 162 def __copy__(self): 163 """Return another optical chain with the same source, optical elements and description-string as this one.""" 164 return OpticalChain(self.source_rays, self.optical_elements, self.detectors, self.description) 165 166 def __deepcopy__(self, memo): 167 """Return another optical chain with the same source, optical elements and description-string as this one.""" 168 return OpticalChain(copy.deepcopy(self.source_rays), copy.deepcopy(self.optical_elements), memo+"\n"+copy.copy(self.description)) 169 170 def get_input_rays(self): 171 """ 172 Returns the list of source rays. 173 """ 174 return self.source_rays 175 176 def get_output_rays(self, force=False, **kwargs): 177 """ 178 Returns the list of (lists of) output rays, calculate them if this hasn't been done yet, 179 or if the source-ray-bundle or anything about the optical elements has changed. 180 181 This is the user-facing method to perform the ray-tracing calculation. 182 """ 183 current_source_rays_hash = hash(self.source_rays) 184 current_optical_elements_hash = mp._hash_list_of_objects(self.optical_elements) 185 if (current_source_rays_hash != self._last_source_rays_hash) or (current_optical_elements_hash != self._last_optical_elements_hash) or force: 186 print("...ray-tracing...", end="", flush=True) 187 self._output_rays = mp.RayTracingCalculation(self.source_rays, self.optical_elements, **kwargs) 188 print( 189 "\r\033[K", end="", flush=True 190 ) # move to beginning of the line with \r and then delete the whole line with \033[K 191 self._last_source_rays_hash = current_source_rays_hash 192 self._last_optical_elements_hash = current_optical_elements_hash 193 194 return self._output_rays 195 196 def get2dPoints(self, Detector = "Focus"): 197 """ 198 Returns the 2D points of the detected rays on the detector. 199 """ 200 if isinstance(Detector, str): 201 if Detector in self.detectors.keys(): 202 Detector = self.detectors[Detector] 203 else: 204 raise ValueError('The "Detector"-argument must be a name of an existing detector or a Detector-object.') 205 return Detector.get_2D_points(self.get_output_rays()[Detector.index])[0] 206 207 def get3dPoints(self, Detector = "Focus"): 208 """ 209 Returns the 3D points of the detected rays on the detector. 210 """ 211 if isinstance(Detector, str): 212 if Detector in self.detectors.keys(): 213 Detector = self.detectors[Detector] 214 else: 215 raise ValueError('The "Detector"-argument must be a name of an existing detector or a Detector-object.') 216 return Detector.get_3D_points(self.get_output_rays()[Detector.index]) 217 218 def getDelays(self, Detector = "Focus"): 219 """ 220 Returns the delays of the detected rays on the detector. 221 """ 222 if isinstance(Detector, str): 223 if Detector in self.detectors.keys(): 224 Detector = self.detectors[Detector] 225 else: 226 raise ValueError('The "Detector"-argument must be a name of an existing detector or a Detector-object.') 227 return Detector.get_Delays(self.get_output_rays()[Detector.index]) 228 229 # %% methods to (mis-)align the optical chain; just uses the corresponding methods of the OpticalElement class... 230 231 def shift_source(self, axis: str|np.ndarray, distance: float): 232 """ 233 Shift source ray bundle by distance (in mm) along the 'axis' specified as 234 a lab-frame vector (numpy-array of length 3) or as one of the strings 235 "vert", "horiz", or "random". 236 237 In the latter case, the reference is the incidence plane of the first 238 non-normal-incidence mirror after the source. If there is none, you will 239 be asked to rather specify the axis as a 3D-numpy-array. 240 241 axis = "vert" means the source position is shifted along the axis perpendicular 242 to that incidence plane, i.e. "vertically" away from the former incidence plane.. 243 244 axis = "horiz" means the source direciton is translated along thr axis in that 245 incidence plane and perpendicular to the current source direction, 246 i.e. "horizontally" in the incidence plane, but retaining the same distance 247 of source and first optical element. 248 249 axis = "random" means the the source direction shifted in a random direction 250 within in the plane perpendicular to the current source direction, 251 e.g. simulating a fluctuation of hte transverse source position. 252 253 Parameters 254 ---------- 255 axis : np.ndarray or str 256 Shift axis, specified either as a 3D lab-frame vector or as one 257 of the strings "vert", "horiz", or "random". 258 259 distance : float 260 Shift distance in mm. 261 262 Returns 263 ------- 264 Nothing, just modifies the property 'source_rays'. 265 """ 266 if type(distance) not in [int, float, np.float64]: 267 raise ValueError('The "distance"-argument must be an int or float number.') 268 269 central_ray_vector = mp.FindCentralRay(self.source_rays).vector 270 mirror_indcs = [i for i, OE in enumerate(self.optical_elements) if "Mirror" in OE.type.type] 271 272 OEnormal = None 273 for i in mirror_indcs: 274 ith_OEnormal = self.optical_elements[i].normal 275 if np.linalg.norm(np.cross(central_ray_vector, ith_OEnormal)) > 1e-10: 276 OEnormal = ith_OEnormal 277 break 278 if OEnormal is None: 279 raise Exception( 280 "There doesn't seem to be a non-normal-incidence mirror in this optical chain, \ 281 so you should rather give 'axis' as a numpy-array of length 3." 282 ) 283 284 if type(axis) == np.ndarray and len(axis) == 3: 285 translation_vector = axis 286 else: 287 perp_axis = np.cross(central_ray_vector, OEnormal) 288 horiz_axis = np.cross(perp_axis, central_ray_vector) 289 290 if axis == "vert": 291 translation_vector = perp_axis 292 elif axis == "horiz": 293 translation_vector = horiz_axis 294 elif axis == "random": 295 translation_vector = ( 296 np.random.uniform(low=-1, high=1, size=1) * perp_axis 297 + np.random.uniform(low=-1, high=1, size=1) * horiz_axis 298 ) 299 else: 300 raise ValueError( 301 'The shift direction must be specified by "axis" as one of ["vert", "horiz", "random"].' 302 ) 303 304 self.source_rays = mgeo.TranslationRayList(self.source_rays, distance * mgeo.Normalize(translation_vector)) 305 306 def tilt_source(self, axis: str|np.ndarray, angle: float): 307 """ 308 Rotate source ray bundle by angle around an axis, specified as 309 a lab-frame vector (numpy-array of length 3) or as one of the strings 310 "in_plane", "out_plane" or "random" direction. 311 312 In the latter case, the function considers the incidence plane of the first 313 non-normal-incidence mirror after the source. If there is none, you will 314 be asked to rather specify the axis as a 3D-numpy-array. 315 316 axis = "in_plane" means the source direction is rotated about an axis 317 perpendicular to that incidence plane, which tilts the source 318 "horizontally" in the same plane. 319 320 axis = "out_plane" means the source direciton is rotated about an axis 321 in that incidence plane and perpendicular to the current source direction, 322 which tilts the source "vertically" out of the former incidence plane. 323 324 axis = "random" means the the source direction is tilted in a random direction, 325 e.g. simulating a beam pointing fluctuation. 326 327 Attention, "angle" is given in deg, so as to remain consitent with the 328 conventions of other functions, although pointing is mostly talked about 329 in mrad instead. 330 331 Parameters 332 ---------- 333 axis : np.ndarray or str 334 Shift axis, specified either as a 3D lab-frame vector or as one 335 of the strings "in_plane", "out_plane", or "random". 336 337 angle : float 338 Rotation angle in degree. 339 340 Returns 341 ------- 342 Nothing, just modifies the property 'source_rays'. 343 """ 344 if type(angle) not in [int, float, np.float64]: 345 raise ValueError('The "angle"-argument must be an int or float number.') 346 347 central_ray_vector = mp.FindCentralRay(self.source_rays).vector 348 mirror_indcs = [i for i, OE in enumerate(self.optical_elements) if "Mirror" in OE.type.type] 349 350 OEnormal = None 351 for i in mirror_indcs: 352 ith_OEnormal = self.optical_elements[i].normal 353 if np.linalg.norm(np.cross(central_ray_vector, ith_OEnormal)) > 1e-10: 354 OEnormal = ith_OEnormal 355 break 356 if OEnormal is None: 357 raise Exception( 358 "There doesn't seem to be a non-normal-incidence mirror in this optical chain, \ 359 so you should rather give 'axis' as a numpy-array of length 3." 360 ) 361 362 if type(axis) == np.ndarray and len(axis) == 3: 363 rot_axis = axis 364 else: 365 rot_axis_in = np.cross(central_ray_vector, OEnormal) 366 rot_axis_out = np.cross(rot_axis_in, central_ray_vector) 367 if axis == "in_plane": 368 rot_axis = rot_axis_in 369 elif axis == "out_plane": 370 rot_axis = rot_axis_out 371 elif axis == "random": 372 rot_axis = ( 373 np.random.uniform(low=-1, high=1, size=1) * rot_axis_in 374 + np.random.uniform(low=-1, high=1, size=1) * rot_axis_out 375 ) 376 else: 377 raise ValueError( 378 'The tilt axis must be specified by as one of ["in_plane", "out_plane", "random"] or as a numpy-array of length 3.' 379 ) 380 381 self.source_rays = mgeo.RotationAroundAxisRayList(self.source_rays, rot_axis, np.deg2rad(angle)) 382 383 def partial_realign(self, OEstart: int, OEstop: int, DistanceList: list, IncidenceAngleList, IncidencePlaneAngleList): 384 """ 385 This as-of-yet not-implemented method will realign only the parts of the optical chain that are between the elements 386 OEstart and OEstop (both included). 387 """ 388 OpticsList = [i.type for i in self.optical_elements[OEstart:OEstop]] 389 outray = self.get_output_rays()[OEstart-1][0] 390 print(outray) 391 new_elements = mp.OEPlacement(OpticsList, DistanceList, IncidenceAngleList, IncidencePlaneAngleList, outray) 392 self.optical_elements[OEstart:OEstop] = new_elements 393 394 # %% 395 def rotate_OE(self, OEindx: int, ref: str, axis: str, angle: float): 396 """ 397 Rotate the optical element OpticalChain.optical_elements[OEindx] with an axis defined relative to the master ray. 398 The axis can either be relative to the incoming master ray or the outgoing one. This is defined by the "ref" variable 399 that takes either of the two values: 400 - "in" 401 - "out" 402 In either case the "axis" can take these values: 403 - "roll": Rotation about the master ray. Looking in the same direction as light propagation, positive is counterclockwise 404 - "pitch": Rotation about the vector in the incidence plane that is normal to the master ray. 405 - "yaw": Rotation about the vector normal to the incidence plane and to the master ray. 406 407 Parameters 408 ---------- 409 OEindx : int 410 Index of the optical element to modify out of OpticalChain.optical_elements. 411 412 ref : str 413 Reference ray used to define the axes of rotation. Can be either: 414 "in" or "out" or "local_normal" (in+out) 415 416 axis : str 417 Rotation axis, specified as one of the strings 418 "pitch", "roll", "yaw" 419 420 angle : float 421 Rotation angle in degree. 422 423 Returns 424 ------- 425 Nothing, just modifies OpticalChain.optical_elements[OEindx]. 426 """ 427 if abs(OEindx) > len(self.optical_elements): 428 raise ValueError( 429 'The "OEnumber"-argument is out of range compared to the length of OpticalChain.optical_elements.' 430 ) 431 if type(angle) not in [int, float, np.float64]: 432 raise ValueError('The "angle"-argument must be an int or float number.') 433 MasterRay = [mp.FindCentralRay(self.source_rays)] 434 TracedRay = mp.RayTracingCalculation(MasterRay, self.optical_elements) 435 TracedRay = [MasterRay] + TracedRay 436 match ref: 437 case "out": 438 RefVec = mgeo.Normalize(TracedRay[OEindx+1][0].vector) 439 case "in": 440 RefVec = mgeo.Normalize(TracedRay[OEindx][0].vector) 441 case "localnormal": 442 In = mgeo.Normalize(TracedRay[OEindx][0].vector) 443 Out = mgeo.Normalize(TracedRay[OEindx+1][0].vector) 444 RefVec = mgeo.Normalize(Out-In) 445 case _: 446 raise ValueError('The "ref"-argument must be a string out of ["in", "out", "localnormal].') 447 if 1 - np.dot(self[OEindx].normal,RefVec) > 1e-2: 448 match axis: 449 case "pitch": 450 self[OEindx].normal = mgeo.RotationAroundAxis(self[OEindx].majoraxis,np.deg2rad(angle),self[OEindx].normal) 451 case "roll": 452 self[OEindx].normal = mgeo.RotationAroundAxis(RefVec,np.deg2rad(angle),self[OEindx].normal) 453 case "yaw": 454 self[OEindx].normal = mgeo.RotationAroundAxis(mgeo.Normalize(np.cross(RefVec, self[OEindx].majoraxis)),np.deg2rad(angle),self[OEindx].normal) 455 case _: 456 raise ValueError('The "axis"-argument must be a string out of ["pitch", "roll", "yaw"].') 457 else: 458 #If the normal vector is aligned with the ray 459 match axis: 460 case "pitch": 461 self[OEindx].normal = mgeo.RotationAroundAxis(self[OEindx].majoraxis,np.deg2rad(angle),self[OEindx].normal) 462 case "roll": 463 self[OEindx].majoraxis = mgeo.RotationAroundAxis(self[OEindx].normal,np.deg2rad(angle),self[OEindx].majoraxis) 464 case "yaw": 465 self[OEindx].normal = mgeo.RotationAroundAxis(mgeo.Normalize(np.cross(RefVec, self[OEindx].majoraxis)),np.deg2rad(angle),self[OEindx].normal) 466 case _: 467 raise ValueError('The "axis"-argument must be a string out of ["pitch", "roll", "yaw"].') 468 469 def shift_OE(self, OEindx: int, ref: str, axis: str, distance: float): 470 """ 471 Shift the optical element OpticalChain.optical_elements[OEindx] along 472 axis referenced to the master ray. 473 474 The axis can either be relative to the incoming master ray or the outgoing one. This is defined by the "ref" variable 475 that takes either of the two values: 476 - "in" 477 - "out" 478 In either case the "axis" can take these values: 479 - "along": Translation along the master ray. 480 - "in_plane": Translation along the vector in the incidence plane that is normal to the master ray. 481 - "out_plane": Translation along the vector normal to the incidence plane and to the master ray. 482 483 Parameters 484 ---------- 485 OEindx : int 486 Index of the optical element to modify out of OpticalChain.optical_elements. 487 488 ref : str 489 Reference ray used to define the axes of rotation. Can be either: 490 "in" or "out". 491 492 axis : str 493 Translation axis, specified as one of the strings 494 "along", "in_plane", "out_plane". 495 496 distance : float 497 Rotation angle in degree. 498 499 Returns 500 ------- 501 Nothing, just modifies OpticalChain.optical_elements[OEindx]. 502 """ 503 if abs(OEindx) >= len(self.optical_elements): 504 raise ValueError( 505 'The "OEnumber"-argument is out of range compared to the length of OpticalChain.optical_elements.' 506 ) 507 if type(distance) not in [int, float, np.float64]: 508 raise ValueError('The "dist"-argument must be an int or float number.') 509 510 MasterRay = [mp.FindCentralRay(self.source_rays)] 511 TracedRay = mp.RayTracingCalculation(MasterRay, self.optical_elements) 512 TracedRay = [MasterRay] + TracedRay 513 match ref: 514 case "out": 515 RefVec = TracedRay[OEindx+1][0].vector 516 case "in": 517 RefVec = TracedRay[OEindx][0].vector 518 case _: 519 raise ValueError('The "ref"-argument must be a string out of ["in", "out"].') 520 match axis: 521 case "along": 522 self[OEindx].position = self[OEindx].position + distance * mgeo.Normalize(RefVec) 523 case "in_plane": 524 self[OEindx].position = self[OEindx].position + distance * mgeo.Normalize(np.cross(RefVec, np.cross(RefVec, self[OEindx].normal))) 525 case "out_plane": 526 self[OEindx].position = self[OEindx].position + distance * mgeo.Normalize(np.cross(RefVec, self[OEindx].normal)) 527 case _: 528 raise ValueError('The "axis"-argument must be a string out of ["along", "in_plane", "out_plane"].') 529 530 # some function that randomly misalings one, or several or all ?
The OpticalChain represents the whole optical setup to be simulated: Its main attributes are a list source-Rays and a list of successive OpticalElements.
The method OpticalChain.get_output_rays() returns an associated list of lists of Rays, each calculated by ray-tracing from one OpticalElement to the next. So OpticalChain.get_output_rays()[i] is the bundle of Rays after optical_elements[i].
The string "description" can contain a short description of the optical setup, or similar notes.
The OpticalChain can be visualized quickly with the method OpticalChain.quickshow(), and more nicely with OpticalChain.render().
The class also provides methods for (mis-)alignment of the source-Ray-bundle and the OpticalElements.
Attributes
source_rays : list[mray.Ray]
List of source rays, which are to be traced.
optical_elements : list[moe.OpticalElement]
List of successive optical elements.
detector: mdet.Detector (optional)
The detector (or list of detectors) to analyse the results.
description : str
A string to describe the optical setup.
Methods
copy_chain()
get_output_rays()
render()
----------
shift_source(axis, distance)
tilt_source(self, axis, angle)
----------
rotate_OE(OEindx, axis, angle)
shift_OE(OEindx, axis, distance)
OpticalChain(source_rays, optical_elements, detectors, description='')
84 def __init__(self, source_rays, optical_elements, detectors, description=""): 85 """ 86 Parameters 87 ---------- 88 source_rays : list[mray.Ray] 89 List of source rays, which are to be traced. 90 91 optical_elements : list[moe.OpticalElement] 92 List of successive optical elements. 93 94 detector: mdet.Detector (optional) 95 The detector (or list of detectors) to analyse the results. 96 97 description : str, optional 98 A string to describe the optical setup. Defaults to ''. 99 100 """ 101 self.source_rays = copy.deepcopy(source_rays) 102 # deepcopy so this object doesn't get changed when the global source_rays change "outside" 103 self.optical_elements = copy.deepcopy(optical_elements) 104 # deepcopy so this object doesn't get changed when the global optical_elements changes "outside" 105 self.detectors = detectors 106 if isinstance(detectors, mdet.Detector): 107 self.detectors = {"Focus": detectors} 108 self.description = description 109 self._output_rays = None 110 self._last_source_rays_hash = None 111 self._last_optical_elements_hash = None
Parameters
source_rays : list[mray.Ray]
List of source rays, which are to be traced.
optical_elements : list[moe.OpticalElement]
List of successive optical elements.
detector: mdet.Detector (optional)
The detector (or list of detectors) to analyse the results.
description : str, optional
A string to describe the optical setup. Defaults to ''.
def
get_input_rays(self):
170 def get_input_rays(self): 171 """ 172 Returns the list of source rays. 173 """ 174 return self.source_rays
Returns the list of source rays.
def
get_output_rays(self, force=False, **kwargs):
176 def get_output_rays(self, force=False, **kwargs): 177 """ 178 Returns the list of (lists of) output rays, calculate them if this hasn't been done yet, 179 or if the source-ray-bundle or anything about the optical elements has changed. 180 181 This is the user-facing method to perform the ray-tracing calculation. 182 """ 183 current_source_rays_hash = hash(self.source_rays) 184 current_optical_elements_hash = mp._hash_list_of_objects(self.optical_elements) 185 if (current_source_rays_hash != self._last_source_rays_hash) or (current_optical_elements_hash != self._last_optical_elements_hash) or force: 186 print("...ray-tracing...", end="", flush=True) 187 self._output_rays = mp.RayTracingCalculation(self.source_rays, self.optical_elements, **kwargs) 188 print( 189 "\r\033[K", end="", flush=True 190 ) # move to beginning of the line with \r and then delete the whole line with \033[K 191 self._last_source_rays_hash = current_source_rays_hash 192 self._last_optical_elements_hash = current_optical_elements_hash 193 194 return self._output_rays
Returns the list of (lists of) output rays, calculate them if this hasn't been done yet, or if the source-ray-bundle or anything about the optical elements has changed.
This is the user-facing method to perform the ray-tracing calculation.
def
get2dPoints(self, Detector='Focus'):
196 def get2dPoints(self, Detector = "Focus"): 197 """ 198 Returns the 2D points of the detected rays on the detector. 199 """ 200 if isinstance(Detector, str): 201 if Detector in self.detectors.keys(): 202 Detector = self.detectors[Detector] 203 else: 204 raise ValueError('The "Detector"-argument must be a name of an existing detector or a Detector-object.') 205 return Detector.get_2D_points(self.get_output_rays()[Detector.index])[0]
Returns the 2D points of the detected rays on the detector.
def
get3dPoints(self, Detector='Focus'):
207 def get3dPoints(self, Detector = "Focus"): 208 """ 209 Returns the 3D points of the detected rays on the detector. 210 """ 211 if isinstance(Detector, str): 212 if Detector in self.detectors.keys(): 213 Detector = self.detectors[Detector] 214 else: 215 raise ValueError('The "Detector"-argument must be a name of an existing detector or a Detector-object.') 216 return Detector.get_3D_points(self.get_output_rays()[Detector.index])
Returns the 3D points of the detected rays on the detector.
def
getDelays(self, Detector='Focus'):
218 def getDelays(self, Detector = "Focus"): 219 """ 220 Returns the delays of the detected rays on the detector. 221 """ 222 if isinstance(Detector, str): 223 if Detector in self.detectors.keys(): 224 Detector = self.detectors[Detector] 225 else: 226 raise ValueError('The "Detector"-argument must be a name of an existing detector or a Detector-object.') 227 return Detector.get_Delays(self.get_output_rays()[Detector.index])
Returns the delays of the detected rays on the detector.
def
shift_source(self, axis: str | numpy.ndarray, distance: float):
231 def shift_source(self, axis: str|np.ndarray, distance: float): 232 """ 233 Shift source ray bundle by distance (in mm) along the 'axis' specified as 234 a lab-frame vector (numpy-array of length 3) or as one of the strings 235 "vert", "horiz", or "random". 236 237 In the latter case, the reference is the incidence plane of the first 238 non-normal-incidence mirror after the source. If there is none, you will 239 be asked to rather specify the axis as a 3D-numpy-array. 240 241 axis = "vert" means the source position is shifted along the axis perpendicular 242 to that incidence plane, i.e. "vertically" away from the former incidence plane.. 243 244 axis = "horiz" means the source direciton is translated along thr axis in that 245 incidence plane and perpendicular to the current source direction, 246 i.e. "horizontally" in the incidence plane, but retaining the same distance 247 of source and first optical element. 248 249 axis = "random" means the the source direction shifted in a random direction 250 within in the plane perpendicular to the current source direction, 251 e.g. simulating a fluctuation of hte transverse source position. 252 253 Parameters 254 ---------- 255 axis : np.ndarray or str 256 Shift axis, specified either as a 3D lab-frame vector or as one 257 of the strings "vert", "horiz", or "random". 258 259 distance : float 260 Shift distance in mm. 261 262 Returns 263 ------- 264 Nothing, just modifies the property 'source_rays'. 265 """ 266 if type(distance) not in [int, float, np.float64]: 267 raise ValueError('The "distance"-argument must be an int or float number.') 268 269 central_ray_vector = mp.FindCentralRay(self.source_rays).vector 270 mirror_indcs = [i for i, OE in enumerate(self.optical_elements) if "Mirror" in OE.type.type] 271 272 OEnormal = None 273 for i in mirror_indcs: 274 ith_OEnormal = self.optical_elements[i].normal 275 if np.linalg.norm(np.cross(central_ray_vector, ith_OEnormal)) > 1e-10: 276 OEnormal = ith_OEnormal 277 break 278 if OEnormal is None: 279 raise Exception( 280 "There doesn't seem to be a non-normal-incidence mirror in this optical chain, \ 281 so you should rather give 'axis' as a numpy-array of length 3." 282 ) 283 284 if type(axis) == np.ndarray and len(axis) == 3: 285 translation_vector = axis 286 else: 287 perp_axis = np.cross(central_ray_vector, OEnormal) 288 horiz_axis = np.cross(perp_axis, central_ray_vector) 289 290 if axis == "vert": 291 translation_vector = perp_axis 292 elif axis == "horiz": 293 translation_vector = horiz_axis 294 elif axis == "random": 295 translation_vector = ( 296 np.random.uniform(low=-1, high=1, size=1) * perp_axis 297 + np.random.uniform(low=-1, high=1, size=1) * horiz_axis 298 ) 299 else: 300 raise ValueError( 301 'The shift direction must be specified by "axis" as one of ["vert", "horiz", "random"].' 302 ) 303 304 self.source_rays = mgeo.TranslationRayList(self.source_rays, distance * mgeo.Normalize(translation_vector))
Shift source ray bundle by distance (in mm) along the 'axis' specified as a lab-frame vector (numpy-array of length 3) or as one of the strings "vert", "horiz", or "random".
In the latter case, the reference is the incidence plane of the first non-normal-incidence mirror after the source. If there is none, you will be asked to rather specify the axis as a 3D-numpy-array.
axis = "vert" means the source position is shifted along the axis perpendicular to that incidence plane, i.e. "vertically" away from the former incidence plane..
axis = "horiz" means the source direciton is translated along thr axis in that incidence plane and perpendicular to the current source direction, i.e. "horizontally" in the incidence plane, but retaining the same distance of source and first optical element.
axis = "random" means the the source direction shifted in a random direction within in the plane perpendicular to the current source direction, e.g. simulating a fluctuation of hte transverse source position.
Parameters
axis : np.ndarray or str
Shift axis, specified either as a 3D lab-frame vector or as one
of the strings "vert", "horiz", or "random".
distance : float
Shift distance in mm.
Returns
Nothing, just modifies the property 'source_rays'.
def
tilt_source(self, axis: str | numpy.ndarray, angle: float):
306 def tilt_source(self, axis: str|np.ndarray, angle: float): 307 """ 308 Rotate source ray bundle by angle around an axis, specified as 309 a lab-frame vector (numpy-array of length 3) or as one of the strings 310 "in_plane", "out_plane" or "random" direction. 311 312 In the latter case, the function considers the incidence plane of the first 313 non-normal-incidence mirror after the source. If there is none, you will 314 be asked to rather specify the axis as a 3D-numpy-array. 315 316 axis = "in_plane" means the source direction is rotated about an axis 317 perpendicular to that incidence plane, which tilts the source 318 "horizontally" in the same plane. 319 320 axis = "out_plane" means the source direciton is rotated about an axis 321 in that incidence plane and perpendicular to the current source direction, 322 which tilts the source "vertically" out of the former incidence plane. 323 324 axis = "random" means the the source direction is tilted in a random direction, 325 e.g. simulating a beam pointing fluctuation. 326 327 Attention, "angle" is given in deg, so as to remain consitent with the 328 conventions of other functions, although pointing is mostly talked about 329 in mrad instead. 330 331 Parameters 332 ---------- 333 axis : np.ndarray or str 334 Shift axis, specified either as a 3D lab-frame vector or as one 335 of the strings "in_plane", "out_plane", or "random". 336 337 angle : float 338 Rotation angle in degree. 339 340 Returns 341 ------- 342 Nothing, just modifies the property 'source_rays'. 343 """ 344 if type(angle) not in [int, float, np.float64]: 345 raise ValueError('The "angle"-argument must be an int or float number.') 346 347 central_ray_vector = mp.FindCentralRay(self.source_rays).vector 348 mirror_indcs = [i for i, OE in enumerate(self.optical_elements) if "Mirror" in OE.type.type] 349 350 OEnormal = None 351 for i in mirror_indcs: 352 ith_OEnormal = self.optical_elements[i].normal 353 if np.linalg.norm(np.cross(central_ray_vector, ith_OEnormal)) > 1e-10: 354 OEnormal = ith_OEnormal 355 break 356 if OEnormal is None: 357 raise Exception( 358 "There doesn't seem to be a non-normal-incidence mirror in this optical chain, \ 359 so you should rather give 'axis' as a numpy-array of length 3." 360 ) 361 362 if type(axis) == np.ndarray and len(axis) == 3: 363 rot_axis = axis 364 else: 365 rot_axis_in = np.cross(central_ray_vector, OEnormal) 366 rot_axis_out = np.cross(rot_axis_in, central_ray_vector) 367 if axis == "in_plane": 368 rot_axis = rot_axis_in 369 elif axis == "out_plane": 370 rot_axis = rot_axis_out 371 elif axis == "random": 372 rot_axis = ( 373 np.random.uniform(low=-1, high=1, size=1) * rot_axis_in 374 + np.random.uniform(low=-1, high=1, size=1) * rot_axis_out 375 ) 376 else: 377 raise ValueError( 378 'The tilt axis must be specified by as one of ["in_plane", "out_plane", "random"] or as a numpy-array of length 3.' 379 ) 380 381 self.source_rays = mgeo.RotationAroundAxisRayList(self.source_rays, rot_axis, np.deg2rad(angle))
Rotate source ray bundle by angle around an axis, specified as a lab-frame vector (numpy-array of length 3) or as one of the strings "in_plane", "out_plane" or "random" direction.
In the latter case, the function considers the incidence plane of the first non-normal-incidence mirror after the source. If there is none, you will be asked to rather specify the axis as a 3D-numpy-array.
axis = "in_plane" means the source direction is rotated about an axis perpendicular to that incidence plane, which tilts the source "horizontally" in the same plane.
axis = "out_plane" means the source direciton is rotated about an axis in that incidence plane and perpendicular to the current source direction, which tilts the source "vertically" out of the former incidence plane.
axis = "random" means the the source direction is tilted in a random direction, e.g. simulating a beam pointing fluctuation.
Attention, "angle" is given in deg, so as to remain consitent with the conventions of other functions, although pointing is mostly talked about in mrad instead.
Parameters
axis : np.ndarray or str
Shift axis, specified either as a 3D lab-frame vector or as one
of the strings "in_plane", "out_plane", or "random".
angle : float
Rotation angle in degree.
Returns
Nothing, just modifies the property 'source_rays'.
def
partial_realign( self, OEstart: int, OEstop: int, DistanceList: list, IncidenceAngleList, IncidencePlaneAngleList):
383 def partial_realign(self, OEstart: int, OEstop: int, DistanceList: list, IncidenceAngleList, IncidencePlaneAngleList): 384 """ 385 This as-of-yet not-implemented method will realign only the parts of the optical chain that are between the elements 386 OEstart and OEstop (both included). 387 """ 388 OpticsList = [i.type for i in self.optical_elements[OEstart:OEstop]] 389 outray = self.get_output_rays()[OEstart-1][0] 390 print(outray) 391 new_elements = mp.OEPlacement(OpticsList, DistanceList, IncidenceAngleList, IncidencePlaneAngleList, outray) 392 self.optical_elements[OEstart:OEstop] = new_elements
This as-of-yet not-implemented method will realign only the parts of the optical chain that are between the elements OEstart and OEstop (both included).
def
rotate_OE(self, OEindx: int, ref: str, axis: str, angle: float):
395 def rotate_OE(self, OEindx: int, ref: str, axis: str, angle: float): 396 """ 397 Rotate the optical element OpticalChain.optical_elements[OEindx] with an axis defined relative to the master ray. 398 The axis can either be relative to the incoming master ray or the outgoing one. This is defined by the "ref" variable 399 that takes either of the two values: 400 - "in" 401 - "out" 402 In either case the "axis" can take these values: 403 - "roll": Rotation about the master ray. Looking in the same direction as light propagation, positive is counterclockwise 404 - "pitch": Rotation about the vector in the incidence plane that is normal to the master ray. 405 - "yaw": Rotation about the vector normal to the incidence plane and to the master ray. 406 407 Parameters 408 ---------- 409 OEindx : int 410 Index of the optical element to modify out of OpticalChain.optical_elements. 411 412 ref : str 413 Reference ray used to define the axes of rotation. Can be either: 414 "in" or "out" or "local_normal" (in+out) 415 416 axis : str 417 Rotation axis, specified as one of the strings 418 "pitch", "roll", "yaw" 419 420 angle : float 421 Rotation angle in degree. 422 423 Returns 424 ------- 425 Nothing, just modifies OpticalChain.optical_elements[OEindx]. 426 """ 427 if abs(OEindx) > len(self.optical_elements): 428 raise ValueError( 429 'The "OEnumber"-argument is out of range compared to the length of OpticalChain.optical_elements.' 430 ) 431 if type(angle) not in [int, float, np.float64]: 432 raise ValueError('The "angle"-argument must be an int or float number.') 433 MasterRay = [mp.FindCentralRay(self.source_rays)] 434 TracedRay = mp.RayTracingCalculation(MasterRay, self.optical_elements) 435 TracedRay = [MasterRay] + TracedRay 436 match ref: 437 case "out": 438 RefVec = mgeo.Normalize(TracedRay[OEindx+1][0].vector) 439 case "in": 440 RefVec = mgeo.Normalize(TracedRay[OEindx][0].vector) 441 case "localnormal": 442 In = mgeo.Normalize(TracedRay[OEindx][0].vector) 443 Out = mgeo.Normalize(TracedRay[OEindx+1][0].vector) 444 RefVec = mgeo.Normalize(Out-In) 445 case _: 446 raise ValueError('The "ref"-argument must be a string out of ["in", "out", "localnormal].') 447 if 1 - np.dot(self[OEindx].normal,RefVec) > 1e-2: 448 match axis: 449 case "pitch": 450 self[OEindx].normal = mgeo.RotationAroundAxis(self[OEindx].majoraxis,np.deg2rad(angle),self[OEindx].normal) 451 case "roll": 452 self[OEindx].normal = mgeo.RotationAroundAxis(RefVec,np.deg2rad(angle),self[OEindx].normal) 453 case "yaw": 454 self[OEindx].normal = mgeo.RotationAroundAxis(mgeo.Normalize(np.cross(RefVec, self[OEindx].majoraxis)),np.deg2rad(angle),self[OEindx].normal) 455 case _: 456 raise ValueError('The "axis"-argument must be a string out of ["pitch", "roll", "yaw"].') 457 else: 458 #If the normal vector is aligned with the ray 459 match axis: 460 case "pitch": 461 self[OEindx].normal = mgeo.RotationAroundAxis(self[OEindx].majoraxis,np.deg2rad(angle),self[OEindx].normal) 462 case "roll": 463 self[OEindx].majoraxis = mgeo.RotationAroundAxis(self[OEindx].normal,np.deg2rad(angle),self[OEindx].majoraxis) 464 case "yaw": 465 self[OEindx].normal = mgeo.RotationAroundAxis(mgeo.Normalize(np.cross(RefVec, self[OEindx].majoraxis)),np.deg2rad(angle),self[OEindx].normal) 466 case _: 467 raise ValueError('The "axis"-argument must be a string out of ["pitch", "roll", "yaw"].')
Rotate the optical element OpticalChain.optical_elements[OEindx] with an axis defined relative to the master ray. The axis can either be relative to the incoming master ray or the outgoing one. This is defined by the "ref" variable that takes either of the two values: - "in" - "out" In either case the "axis" can take these values: - "roll": Rotation about the master ray. Looking in the same direction as light propagation, positive is counterclockwise - "pitch": Rotation about the vector in the incidence plane that is normal to the master ray. - "yaw": Rotation about the vector normal to the incidence plane and to the master ray.
Parameters
OEindx : int
Index of the optical element to modify out of OpticalChain.optical_elements.
ref : str
Reference ray used to define the axes of rotation. Can be either:
"in" or "out" or "local_normal" (in+out)
axis : str
Rotation axis, specified as one of the strings
"pitch", "roll", "yaw"
angle : float
Rotation angle in degree.
Returns
Nothing, just modifies OpticalChain.optical_elements[OEindx].
def
shift_OE(self, OEindx: int, ref: str, axis: str, distance: float):
469 def shift_OE(self, OEindx: int, ref: str, axis: str, distance: float): 470 """ 471 Shift the optical element OpticalChain.optical_elements[OEindx] along 472 axis referenced to the master ray. 473 474 The axis can either be relative to the incoming master ray or the outgoing one. This is defined by the "ref" variable 475 that takes either of the two values: 476 - "in" 477 - "out" 478 In either case the "axis" can take these values: 479 - "along": Translation along the master ray. 480 - "in_plane": Translation along the vector in the incidence plane that is normal to the master ray. 481 - "out_plane": Translation along the vector normal to the incidence plane and to the master ray. 482 483 Parameters 484 ---------- 485 OEindx : int 486 Index of the optical element to modify out of OpticalChain.optical_elements. 487 488 ref : str 489 Reference ray used to define the axes of rotation. Can be either: 490 "in" or "out". 491 492 axis : str 493 Translation axis, specified as one of the strings 494 "along", "in_plane", "out_plane". 495 496 distance : float 497 Rotation angle in degree. 498 499 Returns 500 ------- 501 Nothing, just modifies OpticalChain.optical_elements[OEindx]. 502 """ 503 if abs(OEindx) >= len(self.optical_elements): 504 raise ValueError( 505 'The "OEnumber"-argument is out of range compared to the length of OpticalChain.optical_elements.' 506 ) 507 if type(distance) not in [int, float, np.float64]: 508 raise ValueError('The "dist"-argument must be an int or float number.') 509 510 MasterRay = [mp.FindCentralRay(self.source_rays)] 511 TracedRay = mp.RayTracingCalculation(MasterRay, self.optical_elements) 512 TracedRay = [MasterRay] + TracedRay 513 match ref: 514 case "out": 515 RefVec = TracedRay[OEindx+1][0].vector 516 case "in": 517 RefVec = TracedRay[OEindx][0].vector 518 case _: 519 raise ValueError('The "ref"-argument must be a string out of ["in", "out"].') 520 match axis: 521 case "along": 522 self[OEindx].position = self[OEindx].position + distance * mgeo.Normalize(RefVec) 523 case "in_plane": 524 self[OEindx].position = self[OEindx].position + distance * mgeo.Normalize(np.cross(RefVec, np.cross(RefVec, self[OEindx].normal))) 525 case "out_plane": 526 self[OEindx].position = self[OEindx].position + distance * mgeo.Normalize(np.cross(RefVec, self[OEindx].normal)) 527 case _: 528 raise ValueError('The "axis"-argument must be a string out of ["along", "in_plane", "out_plane"].') 529 530 # some function that randomly misalings one, or several or all ?
Shift the optical element OpticalChain.optical_elements[OEindx] along axis referenced to the master ray.
The axis can either be relative to the incoming master ray or the outgoing one. This is defined by the "ref" variable that takes either of the two values: - "in" - "out" In either case the "axis" can take these values: - "along": Translation along the master ray. - "in_plane": Translation along the vector in the incidence plane that is normal to the master ray. - "out_plane": Translation along the vector normal to the incidence plane and to the master ray.
Parameters
OEindx : int
Index of the optical element to modify out of OpticalChain.optical_elements.
ref : str
Reference ray used to define the axes of rotation. Can be either:
"in" or "out".
axis : str
Translation axis, specified as one of the strings
"along", "in_plane", "out_plane".
distance : float
Rotation angle in degree.
Returns
Nothing, just modifies OpticalChain.optical_elements[OEindx].
def
getETransmission(OpticalChain, IndexIn=0, IndexOut=-1) -> float:
67def _getETransmission(OpticalChain, IndexIn=0, IndexOut=-1) -> float: 68 """ 69 Calculates the energy transmission from the input to the output of the OpticalChain in percent. 70 71 Parameters 72 ---------- 73 OpticalChain : OpticalChain 74 An object of the ModuleOpticalChain.OpticalChain-class. 75 76 IndexIn : int, optional 77 Index of the input RayList in the OpticalChain, defaults to 0. 78 79 IndexOut : int, optional 80 Index of the output RayList in the OpticalChain, defaults to -1. 81 82 Returns 83 ------- 84 ETransmission : float 85 """ 86 Rays = OpticalChain.get_output_rays() 87 if IndexIn == 0: 88 RayListIn = OpticalChain.get_input_rays() 89 else: 90 RayListIn = Rays[IndexIn] 91 RayListOut = Rays[IndexOut] 92 ETransmission = GetETransmission(RayListIn, RayListOut) 93 return ETransmission
Calculates the energy transmission from the input to the output of the OpticalChain in percent.
Parameters
OpticalChain : OpticalChain
An object of the ModuleOpticalChain.OpticalChain-class.
IndexIn : int, optional
Index of the input RayList in the OpticalChain, defaults to 0.
IndexOut : int, optional
Index of the output RayList in the OpticalChain, defaults to -1.
Returns
ETransmission : float
def
getResultsSummary(OpticalChain, Detector='Focus', verbose=False):
147def _getResultsSummary(OpticalChain, Detector = "Focus", verbose=False): 148 """ 149 Calculate and return FocalSpotSize-standard-deviation and Duration-standard-deviation 150 for the given Detector and RayList. 151 If verbose, then also print a summary of the results for the given Detector. 152 153 Parameters 154 ---------- 155 OpticalChain : OpticalChain 156 An object of the ModuleOpticalChain.OpticalChain-class. 157 158 Detector : Detector or str, optional 159 An object of the ModuleDetector.Detector-class or "Focus" to use the focus detector, defaults to "Focus". 160 161 verbose : bool 162 Whether to print a result summary. 163 164 Returns 165 ------- 166 FocalSpotSizeSD : float 167 168 DurationSD : float 169 """ 170 if isinstance(Detector, str): 171 Detector = OpticalChain.detectors[Detector] 172 Index = Detector.index 173 RayListAnalysed = OpticalChain.get_output_rays()[Index] 174 FocalSpotSizeSD, DurationSD = GetResultSummary(Detector, RayListAnalysed, verbose) 175 return FocalSpotSizeSD, DurationSD
Calculate and return FocalSpotSize-standard-deviation and Duration-standard-deviation for the given Detector and RayList. If verbose, then also print a summary of the results for the given Detector.
Parameters
OpticalChain : OpticalChain
An object of the ModuleOpticalChain.OpticalChain-class.
Detector : Detector or str, optional
An object of the ModuleDetector.Detector-class or "Focus" to use the focus detector, defaults to "Focus".
verbose : bool
Whether to print a result summary.
Returns
FocalSpotSizeSD : float
DurationSD : float
def
drawSpotDiagram( OpticalChain, Detector='Focus', DrawAiryAndFourier=False, DrawFocalContour=False, DrawFocal=False, ColorCoded=None, Observer=None) -> matplotlib.figure.Figure:
57def DrawSpotDiagram(OpticalChain, 58 Detector = "Focus", 59 DrawAiryAndFourier=False, 60 DrawFocalContour=False, 61 DrawFocal=False, 62 ColorCoded=None, 63 Observer = None) -> plt.Figure: 64 """ 65 Produce an interactive figure with the spot diagram on the selected Detector. 66 The detector distance can be shifted with the left-right cursor keys. Doing so will actually move the detector. 67 If DrawAiryAndFourier is True, a circle with the Airy-spot-size will be shown. 68 If DrawFocalContour is True, the focal contour calculated from some of the rays will be shown. 69 If DrawFocal is True, a heatmap calculated from some of the rays will be shown. 70 The 'spots' can optionally be color-coded by specifying ColorCoded, which can be one of ["Intensity","Incidence","Delay"]. 71 72 Parameters 73 ---------- 74 RayListAnalysed : list(Ray) 75 List of objects of the ModuleOpticalRay.Ray-class. 76 77 Detector : Detector or str, optional 78 An object of the ModuleDetector.Detector-class or the name of the detector. The default is "Focus". 79 80 DrawAiryAndFourier : bool, optional 81 Whether to draw a circle with the Airy-spot-size. The default is False. 82 83 DrawFocalContour : bool, optional 84 Whether to draw the focal contour. The default is False. 85 86 DrawFocal : bool, optional 87 Whether to draw the focal heatmap. The default is False. 88 89 ColorCoded : str, optional 90 Color-code the spots according to one of ["Intensity","Incidence","Delay"]. The default is None. 91 92 Observer : Observer, optional 93 An observer object. If none, then we just create a copy of the detector and move it when pressing left-right. 94 However, if an observer is specified, then we will change the value of the observer and it will issue 95 the required callbacks to update several plots at the same time. 96 97 Returns 98 ------- 99 fig : matlplotlib-figure-handle. 100 Shows the interactive figure. 101 """ 102 if isinstance(Detector, str): 103 Detector = OpticalChain.detectors[Detector] 104 Index = Detector.index 105 movingDetector = copy(Detector) # We will move this detector when pressing left-right 106 if Observer is None: 107 detectorPosition = Observable(movingDetector.distance) # We will observe the distance of this detector 108 else: 109 detectorPosition = Observer 110 movingDetector.distance = detectorPosition.value 111 112 detectorPosition.register_calculation(lambda x: movingDetector.set_distance(x)) 113 114 RayListAnalysed = OpticalChain.get_output_rays()[Index] 115 116 NumericalAperture = man.GetNumericalAperture(RayListAnalysed, 1) # NA determined from final ray bundle 117 MaxWavelength = np.max([i.wavelength for i in RayListAnalysed]) 118 if DrawAiryAndFourier: 119 AiryRadius = man.GetAiryRadius(MaxWavelength, NumericalAperture) * 1e3 # in µm 120 else: 121 AiryRadius = 0 122 123 if DrawFocalContour or DrawFocal: 124 X,Y,Z = man.GetDiffractionFocus(OpticalChain, movingDetector, Index) 125 Z/=np.max(Z) 126 127 DectectorPoint2D_Xcoord, DectectorPoint2D_Ycoord, FocalSpotSize, SpotSizeSD = mpu._getDetectorPoints( 128 RayListAnalysed, movingDetector 129 ) 130 131 match ColorCoded: 132 case "Intensity": 133 IntensityList = [k.intensity for k in RayListAnalysed] 134 z = np.asarray(IntensityList) 135 zlabel = "Intensity (arb.u.)" 136 title = "Intensity + Spot Diagram\n press left/right to move detector position" 137 addLine = "" 138 case "Incidence": 139 IncidenceList = [np.rad2deg(k.incidence) for k in RayListAnalysed] # degree 140 z = np.asarray(IncidenceList) 141 zlabel = "Incidence angle (deg)" 142 title = "Ray Incidence + Spot Diagram\n press left/right to move detector position" 143 addLine = "" 144 case "Delay": 145 DelayList = movingDetector.get_Delays(RayListAnalysed) 146 IntensityList = [k.intensity for k in RayListAnalysed] 147 DurationSD = mp.WeightedStandardDeviation(DelayList,IntensityList) 148 #DurationSD = mp.StandardDeviation(DelayList) 149 z = np.asarray(DelayList) 150 zlabel = "Delay (fs)" 151 title = "Delay + Spot Diagram\n press left/right to move detector position" 152 addLine = "\n" + "{:.2f}".format(DurationSD) + " fs SD" 153 case _: 154 z = "red" 155 title = "Spot Diagram\n press left/right to move detector position" 156 addLine = "" 157 158 distStep = min(50, max(0.0005, round(FocalSpotSize / 8 / np.arcsin(NumericalAperture) * 10000) / 10000)) # in mm 159 160 plt.ion() 161 fig, ax = plt.subplots() 162 if DrawFocal: 163 focal = ax.pcolormesh(X*1e3,Y*1e3,Z) 164 if DrawFocalContour: 165 levels = [1/np.e**2, 0.5] 166 contour = ax.contourf(X*1e3, Y*1e3, Z, levels=levels, cmap='gray') 167 168 if DrawAiryAndFourier: 169 theta = np.linspace(0, 2 * np.pi, 100) 170 x = AiryRadius * np.cos(theta) 171 y = AiryRadius * np.sin(theta) # 172 ax.plot(x, y, c="black") 173 174 175 foo = ax.scatter( 176 DectectorPoint2D_Xcoord, 177 DectectorPoint2D_Ycoord, 178 c=z, 179 s=15, 180 label="{:.3f}".format(detectorPosition.value) + " mm\n" + "{:.1f}".format(SpotSizeSD * 1e3) + " \u03BCm SD" + addLine, 181 ) 182 183 axisLim = 1.1 * max(AiryRadius, 0.5 * FocalSpotSize * 1000) 184 ax.set_xlim(-axisLim, axisLim) 185 ax.set_ylim(-axisLim, axisLim) 186 187 if ColorCoded == "Intensity" or ColorCoded == "Incidence" or ColorCoded == "Delay": 188 cbar = fig.colorbar(foo) 189 cbar.set_label(zlabel) 190 191 ax.legend(loc="upper right") 192 ax.set_xlabel("X (µm)") 193 ax.set_ylabel("Y (µm)") 194 ax.set_title(title) 195 # ax.margins(x=0) 196 197 198 def update_plot(new_value): 199 nonlocal movingDetector, ColorCoded, zlabel, cbar, detectorPosition, foo, distStep, focal, contour, levels, Index, RayListAnalysed 200 201 newDectectorPoint2D_Xcoord, newDectectorPoint2D_Ycoord, newFocalSpotSize, newSpotSizeSD = mpu._getDetectorPoints( 202 RayListAnalysed, movingDetector 203 ) 204 205 if DrawFocal: 206 focal.set_array(Z) 207 if DrawFocalContour: 208 levels = [1/np.e**2, 0.5] 209 for coll in contour.collections: 210 coll.remove() # Remove old contour lines 211 contour = ax.contourf(X * 1e3, Y * 1e3, Z, levels=levels, cmap='gray') 212 213 xy = foo.get_offsets() 214 xy[:, 0] = newDectectorPoint2D_Xcoord 215 xy[:, 1] = newDectectorPoint2D_Ycoord 216 foo.set_offsets(xy) 217 218 219 if ColorCoded == "Delay": 220 newDelayList = np.asarray(movingDetector.get_Delays(RayListAnalysed)) 221 newDurationSD = mp.WeightedStandardDeviation(newDelayList,IntensityList) 222 #newDurationSD = mp.StandardDeviation(newDelayList) 223 newaddLine = "\n" + "{:.2f}".format(newDurationSD) + " fs SD" 224 foo.set_array(newDelayList) 225 foo.set_clim(min(newDelayList), max(newDelayList)) 226 cbar.update_normal(foo) 227 else: 228 newaddLine = "" 229 230 foo.set_label( 231 "{:.3f}".format(detectorPosition.value) + " mm\n" + "{:.1f}".format(newSpotSizeSD * 1e3) + " \u03BCm SD" + newaddLine 232 ) 233 ax.legend(loc="upper right") 234 235 axisLim = 1.1 * max(AiryRadius, 0.5 * newFocalSpotSize * 1000) 236 ax.set_xlim(-axisLim, axisLim) 237 ax.set_ylim(-axisLim, axisLim) 238 239 distStep = min( 240 50, max(0.0005, round(newFocalSpotSize / 8 / np.arcsin(NumericalAperture) * 10000) / 10000) 241 ) # in mm 242 243 fig.canvas.draw_idle() 244 245 246 def press(event): 247 nonlocal detectorPosition, distStep 248 if event.key == "right": 249 detectorPosition.value += distStep 250 elif event.key == "left": 251 if detectorPosition.value > 1.5 * distStep: 252 detectorPosition.value -= distStep 253 else: 254 detectorPosition.value = 0.5 * distStep 255 else: 256 return None 257 258 fig.canvas.mpl_connect("key_press_event", press) 259 260 plt.show() 261 262 detectorPosition.register(update_plot) 263 264 265 return fig, detectorPosition
Produce an interactive figure with the spot diagram on the selected Detector. The detector distance can be shifted with the left-right cursor keys. Doing so will actually move the detector. If DrawAiryAndFourier is True, a circle with the Airy-spot-size will be shown. If DrawFocalContour is True, the focal contour calculated from some of the rays will be shown. If DrawFocal is True, a heatmap calculated from some of the rays will be shown. The 'spots' can optionally be color-coded by specifying ColorCoded, which can be one of ["Intensity","Incidence","Delay"].
Parameters
RayListAnalysed : list(Ray)
List of objects of the ModuleOpticalRay.Ray-class.
Detector : Detector or str, optional
An object of the ModuleDetector.Detector-class or the name of the detector. The default is "Focus".
DrawAiryAndFourier : bool, optional
Whether to draw a circle with the Airy-spot-size. The default is False.
DrawFocalContour : bool, optional
Whether to draw the focal contour. The default is False.
DrawFocal : bool, optional
Whether to draw the focal heatmap. The default is False.
ColorCoded : str, optional
Color-code the spots according to one of ["Intensity","Incidence","Delay"]. The default is None.
Observer : Observer, optional
An observer object. If none, then we just create a copy of the detector and move it when pressing left-right.
However, if an observer is specified, then we will change the value of the observer and it will issue
the required callbacks to update several plots at the same time.
Returns
fig : matlplotlib-figure-handle.
Shows the interactive figure.
def
drawDelaySpots( OpticalChain, DeltaFT: tuple[int, float], Detector='Focus', DrawAiryAndFourier=False, ColorCoded=None, Observer=None) -> matplotlib.figure.Figure:
270def DrawDelaySpots(OpticalChain, 271 DeltaFT: tuple[int, float], 272 Detector = "Focus", 273 DrawAiryAndFourier=False, 274 ColorCoded=None, 275 Observer = None 276 ) -> plt.Figure: 277 """ 278 Produce a an interactive figure with a spot diagram resulting from the RayListAnalysed 279 hitting the Detector, with the ray-delays shown in the 3rd dimension. 280 The detector distance can be shifted with the left-right cursor keys. 281 If DrawAiryAndFourier is True, a cylinder is shown whose diameter is the Airy-spot-size and 282 whose height is the Fourier-limited pulse duration 'given by 'DeltaFT'. 283 284 The 'spots' can optionally be color-coded by specifying ColorCoded as ["Intensity","Incidence"]. 285 286 Parameters 287 ---------- 288 RayListAnalysed : list(Ray) 289 List of objects of the ModuleOpticalRay.Ray-class. 290 291 Detector : Detector 292 An object of the ModuleDetector.Detector-class. 293 294 DeltaFT : (int, float) 295 The Fourier-limited pulse duration. Just used as a reference to compare the temporal spread 296 induced by the ray-delays. 297 298 DrawAiryAndFourier : bool, optional 299 Whether to draw a cylinder showing the Airy-spot-size and Fourier-limited-duration. 300 The default is False. 301 302 ColorCoded : str, optional 303 Color-code the spots according to one of ["Intensity","Incidence"]. 304 The default is None. 305 306 Returns 307 ------- 308 fig : matlplotlib-figure-handle. 309 Shows the interactive figure. 310 """ 311 if isinstance(Detector, str): 312 Det = OpticalChain.detectors[Detector] 313 else: 314 Det = Detector 315 Index = Det.index 316 Detector = copy(Det) 317 if Observer is None: 318 detectorPosition = Observable(Detector.distance) 319 else: 320 detectorPosition = Observer 321 Detector.distance = detectorPosition.value 322 323 RayListAnalysed = OpticalChain.get_output_rays()[Index] 324 fig, NumericalAperture, AiryRadius, FocalSpotSize = _drawDelayGraph( 325 RayListAnalysed, Detector, detectorPosition.value, DeltaFT, DrawAiryAndFourier, ColorCoded 326 ) 327 328 distStep = min(50, max(0.0005, round(FocalSpotSize / 8 / np.arcsin(NumericalAperture) * 10000) / 10000)) # in mm 329 330 movingDetector = copy(Detector) 331 332 def update_plot(new_value): 333 nonlocal movingDetector, ColorCoded, detectorPosition, distStep, fig 334 ax = fig.axes[0] 335 cam = [ax.azim, ax.elev, ax._dist] 336 fig, sameNumericalAperture, sameAiryRadius, newFocalSpotSize = _drawDelayGraph( 337 RayListAnalysed, movingDetector, detectorPosition.value, DeltaFT, DrawAiryAndFourier, ColorCoded, fig 338 ) 339 ax = fig.axes[0] 340 ax.azim, ax.elev, ax._dist = cam 341 distStep = min( 342 50, max(0.0005, round(newFocalSpotSize / 8 / np.arcsin(NumericalAperture) * 10000) / 10000) 343 ) 344 return fig 345 346 def press(event): 347 nonlocal detectorPosition, distStep, movingDetector, fig 348 if event.key == "right": 349 detectorPosition.value += distStep 350 elif event.key == "left": 351 if detectorPosition.value > 1.5 * distStep: 352 detectorPosition.value -= distStep 353 else: 354 detectorPosition.value = 0.5 * distStep 355 356 fig.canvas.mpl_connect("key_press_event", press) 357 detectorPosition.register(update_plot) 358 detectorPosition.register_calculation(lambda x: movingDetector.set_distance(x)) 359 360 return fig, Observable
Produce a an interactive figure with a spot diagram resulting from the RayListAnalysed hitting the Detector, with the ray-delays shown in the 3rd dimension. The detector distance can be shifted with the left-right cursor keys. If DrawAiryAndFourier is True, a cylinder is shown whose diameter is the Airy-spot-size and whose height is the Fourier-limited pulse duration 'given by 'DeltaFT'.
The 'spots' can optionally be color-coded by specifying ColorCoded as ["Intensity","Incidence"].
Parameters
RayListAnalysed : list(Ray)
List of objects of the ModuleOpticalRay.Ray-class.
Detector : Detector
An object of the ModuleDetector.Detector-class.
DeltaFT : (int, float)
The Fourier-limited pulse duration. Just used as a reference to compare the temporal spread
induced by the ray-delays.
DrawAiryAndFourier : bool, optional
Whether to draw a cylinder showing the Airy-spot-size and Fourier-limited-duration.
The default is False.
ColorCoded : str, optional
Color-code the spots according to one of ["Intensity","Incidence"].
The default is None.
Returns
fig : matlplotlib-figure-handle.
Shows the interactive figure.
def
drawMirrorProjection( OpticalChain, ReflectionNumber: int, ColorCoded=None, Detector='') -> matplotlib.figure.Figure:
441def DrawMirrorProjection(OpticalChain, ReflectionNumber: int, ColorCoded=None, Detector="") -> plt.Figure: 442 """ 443 Produce a plot of the ray impact points on the optical element with index 'ReflectionNumber'. 444 The points can be color-coded according ["Incidence","Intensity","Delay"], where the ray delay is 445 measured at the Detector. 446 447 Parameters 448 ---------- 449 OpticalChain : OpticalChain 450 List of objects of the ModuleOpticalOpticalChain.OpticalChain-class. 451 452 ReflectionNumber : int 453 Index specifying the optical element on which you want to see the impact points. 454 455 Detector : Detector, optional 456 Object of the ModuleDetector.Detector-class. Only necessary to project delays. The default is None. 457 458 ColorCoded : str, optional 459 Specifies which ray property to color-code: ["Incidence","Intensity","Delay"]. The default is None. 460 461 Returns 462 ------- 463 fig : matlplotlib-figure-handle. 464 Shows the figure. 465 """ 466 from mpl_toolkits.axes_grid1 import make_axes_locatable 467 if isinstance(Detector, str): 468 if Detector == "": 469 Detector = None 470 else: 471 Detector = OpticalChain.detectors[Detector] 472 473 Position = OpticalChain[ReflectionNumber].position 474 q = OpticalChain[ReflectionNumber].orientation 475 # n = OpticalChain.optical_elements[ReflectionNumber].normal 476 # m = OpticalChain.optical_elements[ReflectionNumber].majoraxis 477 478 RayListAnalysed = OpticalChain.get_output_rays()[ReflectionNumber] 479 # transform rays into the mirror-support reference frame 480 # (same as mirror frame but without the shift by mirror-centre) 481 r0 = OpticalChain[ReflectionNumber].r0 482 RayList = [r.to_basis(*OpticalChain[ReflectionNumber].basis) for r in RayListAnalysed] 483 484 x = np.asarray([k.point[0] for k in RayList]) - r0[0] 485 y = np.asarray([k.point[1] for k in RayList]) - r0[1] 486 if ColorCoded == "Intensity": 487 IntensityList = [k.intensity for k in RayListAnalysed] 488 z = np.asarray(IntensityList) 489 zlabel = "Intensity (arb.u.)" 490 title = "Ray intensity projected on mirror " 491 elif ColorCoded == "Incidence": 492 IncidenceList = [np.rad2deg(k.incidence) for k in RayListAnalysed] # in degree 493 z = np.asarray(IncidenceList) 494 zlabel = "Incidence angle (deg)" 495 title = "Ray incidence projected on mirror " 496 elif ColorCoded == "Delay": 497 if Detector is not None: 498 z = np.asarray(Detector.get_Delays(RayListAnalysed)) 499 zlabel = "Delay (fs)" 500 title = "Ray delay at detector projected on mirror " 501 else: 502 raise ValueError("If you want to project ray delays, you must specify a detector.") 503 else: 504 z = "red" 505 title = "Ray impact points projected on mirror" 506 507 plt.ion() 508 fig = plt.figure() 509 ax = OpticalChain.optical_elements[ReflectionNumber].support._ContourSupport(fig) 510 p = plt.scatter(x, y, c=z, s=15) 511 if ColorCoded == "Delay" or ColorCoded == "Incidence" or ColorCoded == "Intensity": 512 divider = make_axes_locatable(ax) 513 cax = divider.append_axes("right", size="5%", pad=0.05) 514 cbar = fig.colorbar(p, cax=cax) 515 cbar.set_label(zlabel) 516 ax.set_xlabel("x (mm)") 517 ax.set_ylabel("y (mm)") 518 plt.title(title, loc="right") 519 plt.tight_layout() 520 521 bbox = ax.get_position() 522 bbox.set_points(bbox.get_points() - np.array([[0.01, 0], [0.01, 0]])) 523 ax.set_position(bbox) 524 plt.show() 525 526 return fig
Produce a plot of the ray impact points on the optical element with index 'ReflectionNumber'. The points can be color-coded according ["Incidence","Intensity","Delay"], where the ray delay is measured at the Detector.
Parameters
OpticalChain : OpticalChain
List of objects of the ModuleOpticalOpticalChain.OpticalChain-class.
ReflectionNumber : int
Index specifying the optical element on which you want to see the impact points.
Detector : Detector, optional
Object of the ModuleDetector.Detector-class. Only necessary to project delays. The default is None.
ColorCoded : str, optional
Specifies which ray property to color-code: ["Incidence","Intensity","Delay"]. The default is None.
Returns
fig : matlplotlib-figure-handle.
Shows the figure.
def
render( OpticalChain, EndDistance=None, maxRays=300, OEpoints=2000, draw_mesh=False, cycle_ray_colors=False, impact_points=False, DrawDetectors=True, DetectedRays=False, Observers={}):
532def DrawSetup(OpticalChain, 533 EndDistance=None, 534 maxRays=300, 535 OEpoints=2000, 536 draw_mesh=False, 537 cycle_ray_colors = False, 538 impact_points = False, 539 DrawDetectors=True, 540 DetectedRays = False, 541 Observers = dict()): 542 """ 543 Renders an image of the Optical setup and the traced rays. 544 545 Parameters 546 ---------- 547 OpticalChain : OpticalChain 548 List of objects of the ModuleOpticalOpticalChain.OpticalChain-class. 549 550 EndDistance : float, optional 551 The rays of the last ray bundle are drawn with a length given by EndDistance (in mm). If not specified, 552 this distance is set to that between the source point and the 1st optical element. 553 554 maxRays: int 555 The maximum number of rays to render. Rendering all the traced rays is a insufferable resource hog 556 and not required for a nice image. Default is 150. 557 558 OEpoints : int 559 How many little spheres to draw to represent the optical elements. Default is 2000. 560 561 Returns 562 ------- 563 fig : Pyvista-figure-handle. 564 Shows the figure. 565 """ 566 567 RayListHistory = [OpticalChain.source_rays] + OpticalChain.get_output_rays() 568 569 if EndDistance is None: 570 EndDistance = np.linalg.norm(OpticalChain.source_rays[0].point - OpticalChain.optical_elements[0].position) 571 572 print("...rendering image of optical chain...", end="", flush=True) 573 fig = pvqt.BackgroundPlotter(window_size=(1500, 500), notebook=False) # Opening a window 574 fig.set_background('white') 575 576 if cycle_ray_colors: 577 colors = mpu.generate_distinct_colors(len(OpticalChain)+1) 578 else: 579 colors = [[0.7, 0, 0]]*(len(OpticalChain)+1) # Default color: dark red 580 581 # Optics display 582 # For each optic we will send the figure to the function _RenderOpticalElement and it will add the optic to the figure 583 for i,OE in enumerate(OpticalChain.optical_elements): 584 color = pv.Color(colors[i+1]) 585 rgb = color.float_rgb 586 h, l, s = rgb_to_hls(*rgb) 587 s = max(0, min(1, s * 0.3)) # Decrease saturation 588 l = max(0, min(1, l + 0.1)) # Increase lightness 589 new_rgb = hls_to_rgb(h, l, s) 590 darkened_color = pv.Color(new_rgb) 591 mpm._RenderOpticalElement(fig, OE, OEpoints, draw_mesh, darkened_color, index=i) 592 ray_meshes = mpm._RenderRays(RayListHistory, EndDistance, maxRays) 593 for i,ray in enumerate(ray_meshes): 594 color = pv.Color(colors[i]) 595 fig.add_mesh(ray, color=color, name=f"RayBundle_{i}") 596 if impact_points: 597 for i,rays in enumerate(RayListHistory): 598 points = np.array([list(r.point) for r in rays], dtype=np.float32) 599 points = pv.PolyData(points) 600 color = pv.Color(colors[i-1]) 601 fig.add_mesh(points, color=color, point_size=5, name=f"RayImpactPoints_{i}") 602 603 detector_copies = {key: copy(OpticalChain.detectors[key]) for key in OpticalChain.detectors.keys()} 604 detector_meshes_list = [] 605 detectedpoint_meshes = dict() 606 607 if OpticalChain.detectors is not None and DrawDetectors: 608 # Detector display 609 for key in OpticalChain.detectors.keys(): 610 det = detector_copies[key] 611 index = OpticalChain.detectors[key].index 612 if key in Observers: 613 det.distance = Observers[key].value 614 #Observers[key].register_calculation(lambda x: det.set_distance(x)) 615 mpm._RenderDetector(fig, det, name = key, detector_meshes = detector_meshes_list) 616 if DetectedRays: 617 RayListAnalysed = OpticalChain.get_output_rays()[index] 618 points = det.get_3D_points(RayListAnalysed) 619 points = pv.PolyData(points) 620 detectedpoint_meshes[key] = points 621 fig.add_mesh(points, color='purple', point_size=5, name=f"DetectedRays_{key}") 622 detector_meshes = dict(zip(OpticalChain.detectors.keys(), detector_meshes_list)) 623 624 # Now we define a function that will move on the plot the detector with name "detname" when it's called 625 def move_detector(detname, new_value): 626 nonlocal fig, detector_meshes, detectedpoint_meshes, DetectedRays, detectedpoint_meshes, detector_copies, OpticalChain 627 det = detector_copies[detname] 628 index = OpticalChain.detectors[detname].index 629 det_mesh = detector_meshes[detname] 630 translation = det.normal * (det.distance - new_value) 631 det_mesh.translate(translation, inplace=True) 632 det.distance = new_value 633 if DetectedRays: 634 points_mesh = detectedpoint_meshes[detname] 635 points_mesh.points = det.get_3D_points(OpticalChain.get_output_rays()[index]) 636 fig.show() 637 638 # Now we register the function to the observers 639 for key in OpticalChain.detectors.keys(): 640 if key in Observers: 641 Observers[key].register(lambda x: move_detector(key, x)) 642 643 #pv.save_meshio('optics.inp', pointcloud) 644 print( 645 "\r\033[K", end="", flush=True 646 ) # move to beginning of the line with \r and then delete the whole line with \033[K 647 fig.show() 648 return fig
Renders an image of the Optical setup and the traced rays.
Parameters
OpticalChain : OpticalChain
List of objects of the ModuleOpticalOpticalChain.OpticalChain-class.
EndDistance : float, optional
The rays of the last ray bundle are drawn with a length given by EndDistance (in mm). If not specified,
this distance is set to that between the source point and the 1st optical element.
maxRays: int
The maximum number of rays to render. Rendering all the traced rays is a insufferable resource hog
and not required for a nice image. Default is 150.
OEpoints : int
How many little spheres to draw to represent the optical elements. Default is 2000.
Returns
fig : Pyvista-figure-handle.
Shows the figure.
def
drawCaustics(OpticalChain, Range=1, Detector='Focus', Npoints=1000, Nrays=1000):
704def DrawCaustics(OpticalChain, Range=1, Detector="Focus" , Npoints=1000, Nrays=1000): 705 """ 706 This function displays the caustics of the rays on the detector. 707 To do so, it calculates the intersections of the rays with the detector over a 708 range determined by the parameter Range, and then plots the standard deviation of the 709 positions in the x and y directions. 710 711 Parameters 712 ---------- 713 OpticalChain : OpticalChain 714 The optical chain to analyse. 715 716 DetectorName : str 717 The name of the detector on which the caustics are calculated. 718 719 Range : float 720 The range of the detector over which to calculate the caustics. 721 722 Npoints : int, optional 723 The number of points to sample on the detector. The default is 1000. 724 725 Returns 726 ------- 727 fig : Figure 728 The figure of the plot. 729 """ 730 distances = np.linspace(-Range, Range, Npoints) 731 if isinstance(Detector, str): 732 Det = OpticalChain.detectors[Detector] 733 Index = Det.index 734 Rays = OpticalChain.get_output_rays()[Index] 735 Nrays = min(Nrays, len(Rays)) 736 Rays = np.random.choice(Rays, Nrays, replace=False) 737 LocalRayList = [r.to_basis(*Det.basis) for r in Rays] 738 Points = mgeo.IntersectionRayListZPlane(LocalRayList, distances) 739 x_std = [] 740 y_std = [] 741 for i in range(len(distances)): 742 x_std.append(mp.StandardDeviation(Points[i][:,0])) 743 y_std.append(mp.StandardDeviation(Points[i][:,1])) 744 plt.ion() 745 fig, ax = plt.subplots() 746 ax.plot(distances, x_std, label="x std") 747 ax.plot(distances, y_std, label="y std") 748 ax.set_xlabel("Detector distance (mm)") 749 ax.set_ylabel("Standard deviation (mm)") 750 ax.legend() 751 plt.title("Caustics") 752 plt.show() 753 return fig
This function displays the caustics of the rays on the detector. To do so, it calculates the intersections of the rays with the detector over a range determined by the parameter Range, and then plots the standard deviation of the positions in the x and y directions.
Parameters
OpticalChain : OpticalChain The optical chain to analyse.
DetectorName : str The name of the detector on which the caustics are calculated.
Range : float The range of the detector over which to calculate the caustics.
Npoints : int, optional The number of points to sample on the detector. The default is 1000.
Returns
fig : Figure The figure of the plot.