|
def | tel_init (carmaWrap_context context, conf.Param_geom p_geom, conf.Param_tel p_tel, r0=None, ittime=None, p_wfss=None, dm=None) |
| Initialize the overall geometry of the AO system, including pupil and WFS. More...
|
|
def | init_wfs_geom (conf.Param_wfs p_wfs, float r0, conf.Param_tel p_tel, conf.Param_geom p_geom, float ittime, verbose=1) |
| Compute the geometry of WFSs: valid subaps, positions of the subaps, flux per subap, etc... More...
|
|
def | init_wfs_size (conf.Param_wfs p_wfs, float r0, conf.Param_tel p_tel, verbose=1) |
| Compute all the parameters usefull for further WFS image computation (array sizes) More...
|
|
def | compute_nphotons (wfs_type, ittime, optthroughput, diam, cobs=0, nxsub=0, zerop=0, gsmag=0, lgsreturnperwatt=0, laserpower=0, verbose=1) |
| Determines the number of photons TBC. More...
|
|
def | init_pyrhr_geom (conf.Param_wfs p_wfs, float r0, conf.Param_tel p_tel, conf.Param_geom p_geom, float ittime, bool verbose=True) |
| Compute the geometry of PYRHR WFSs: valid subaps, positions of the subaps, flux per subap, etc... More...
|
|
def | init_sh_geom (conf.Param_wfs p_wfs, float r0, conf.Param_tel p_tel, conf.Param_geom p_geom, float ittime, bool verbose=True) |
| Compute the geometry of SH WFSs: valid subaps, positions of the subaps, flux per subap, etc... More...
|
|
def | geom_init (conf.Param_geom p_geom, conf.Param_tel p_tel, padding=2) |
| Initialize the system geometry. More...
|
|
def | geom_init_generic (p_geom, pupdiam, t_spiders=0.01, spiders_type="six", xc=0, yc=0, real=0, cobs=0) |
| Initialize the system geometry. More...
|
|
def | pad_array (A, N) |
|
Initialization of the system geometry and of the Telescope object.
- Author
- COMPASS Team https://github.com/ANR-COMPASS
- Version
- 5.0.0
- Date
- 2020/05/18
- Copyright
- GNU Lesser General Public License
This file is part of COMPASS https://anr-compass.github.io/compass/
Copyright (C) 2011-2019 COMPASS Team https://github.com/ANR-COMPASS All rights reserved. Distributed under GNU - LGPL
COMPASS is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or any later version.
COMPASS: End-to-end AO simulation tool using GPU acceleration The COMPASS platform was designed to meet the need of high-performance for the simulation of AO systems.
The final product includes a software package for simulating all the critical subcomponents of AO, particularly in the context of the ELT and a real-time core based on several control approaches, with performances consistent with its integration into an instrument. Taking advantage of the specific hardware architecture of the GPU, the COMPASS tool allows to achieve adequate execution speeds to conduct large simulation campaigns called to the ELT.
The COMPASS platform can be used to carry a wide variety of simulations to both testspecific components of AO of the E-ELT (such as wavefront analysis device with a pyramid or elongated Laser star), and various systems configurations such as multi-conjugate AO.
COMPASS is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public License along with COMPASS. If not, see https://www.gnu.org/licenses/lgpl-3.0.txt.
def shesha.init.geom_init.compute_nphotons |
( |
|
wfs_type, |
|
|
|
ittime, |
|
|
|
optthroughput, |
|
|
|
diam, |
|
|
|
cobs = 0 , |
|
|
|
nxsub = 0 , |
|
|
|
zerop = 0 , |
|
|
|
gsmag = 0 , |
|
|
|
lgsreturnperwatt = 0 , |
|
|
|
laserpower = 0 , |
|
|
|
verbose = 1 |
|
) |
| |
Determines the number of photons TBC.
:parameters: wfs_type (scons.WFSType) : wfs type: SH or PYRHR.
ittime (float) : 1/loop frequency [s].
optthroughput (float) : wfs global throughput.
diam (float) : telescope diameter.
cobs (float) : (optional for SH) telescope central obstruction.
nxsub (int) : (optional for PYRHR) linear number of subaps.
zerop (float) : (optional for LGS) detector zero point expressed in ph/m**2/s in the bandwidth of the WFS.
gsmag (float) : (optional for LGS) magnitude of guide star.
lgsreturnperwatt (float) : (optional for NGS) return per watt factor (high season : 10 ph/cm2/s/W).
laserpower (float) : (optional for NGS) laser power in W.
verbose (bool) : (optional) display informations if True.
for PYRHR WFS: nphotons = compute_nphotons(scons.WFSType.PYRHR, ittime, optthroughput, diam, cobs=?, zerop=?, gsmag=?) for NGS SH WFS: nphotons = compute_nphotons(scons.WFSType.SH, ittime, optthroughput, diam, nxsub=?, zerop=?, gsmag=?) for LGS SH WFS: nphotons = compute_nphotons(scons.WFSType.SH, ittime, optthroughput, diam, nxsub=?, lgsreturnperwatt=?, laserpower=?)
Definition at line 376 of file geom_init.py.
376 lgsreturnperwatt=?, laserpower=?)
380 if (wfs_type == scons.WFSType.PYRHR or wfs_type == scons.WFSType.PYRLR):
381 surface = np.pi / 4. * (1 - cobs**2.) * diam**2.
382 elif (wfs_type == scons.WFSType.SH):
383 # from the guide star
384 if (laserpower == 0):
387 surface = (diam / nxsub)**2.
388 # include throughput to WFS for unobstructed
389 # subaperture per iteration
390 else: # we are dealing with a LGS
391 nphotons = lgsreturnperwatt * laserpower * \
392 optthroughput * (diam / nxsub) ** 2. * 1e4 * ittime
394 # ... for given power include throughput to WFS
395 # for unobstructed subaperture per iteration
397 print("nphotons : ", nphotons)
400 raise RuntimeError("WFS unknown")
402 nphotons = zerop * 10. ** (-0.4 * gsmag) * ittime * \
403 optthroughput * surface
406 print("nphotons : ", nphotons)
def shesha.init.geom_init.init_pyrhr_geom |
( |
conf.Param_wfs |
p_wfs, |
|
|
float |
r0, |
|
|
conf.Param_tel |
p_tel, |
|
|
conf.Param_geom |
p_geom, |
|
|
float |
ittime, |
|
|
bool |
verbose = True |
|
) |
| |
Compute the geometry of PYRHR WFSs: valid subaps, positions of the subaps, flux per subap, etc...
:parameters: p_wfs (Param_wfs) : wfs settings
r0 (float) : atmos r0 @ 0.5 microns
p_tel (Param_tel) : telescope settings
geom (Param_geom) : geom settings
ittime (float) : 1/loop frequency [s]
verbose (bool) : (optional) display informations if True
Definition at line 427 of file geom_init.py.
430 p_wfs.npix = p_wfs._pdiam
431 p_wfs._Ntot = p_geom._n
433 # Creating field stop mask
434 fsradius_pixels = int(p_wfs.fssize / p_wfs._qpixsize / 2.)
435 if (p_wfs.fstop == scons.FieldStopType.ROUND):
436 focmask = util.dist(p_wfs._Nfft, xc=p_wfs._Nfft / 2. - 0.5,
437 yc=p_wfs._Nfft / 2. - 0.5) < (fsradius_pixels)
438 elif (p_wfs.fstop == scons.FieldStopType.SQUARE):
439 X = np.indices((p_wfs._Nfft, p_wfs._Nfft)) + 1 # TODO: +1 ??
440 x = X[1] - (p_wfs._Nfft + 1.) / 2.
441 y = X[0] - (p_wfs._Nfft + 1.) / 2.
442 focmask = (np.abs(x) <= (fsradius_pixels)) * \
443 (np.abs(y) <= (fsradius_pixels))
445 msg = "PYRHR wfs fstop must be FieldStopType.[ROUND|SQUARE]"
446 raise ValueError(msg)
448 pyr_focmask = focmask * 1.0 # np.fft.fftshift(focmask*1.0)
449 p_wfs._submask = np.fft.fftshift(pyr_focmask)
451 # Creating pyramid mask
452 pyrsize = p_wfs._Nfft
454 rpup = p_geom.pupdiam / 2.0
455 dpup = p_geom.pupdiam
456 nrebin = p_wfs._nrebin
457 fracsub = p_wfs.fracsub
458 if p_wfs.pyr_pup_sep == -1:
459 pup_sep = int(p_wfs.nxsub)
461 pup_sep = p_wfs.pyr_pup_sep
462 # number of pix between two centers two pupil images
464 y = np.tile(np.arange(pyrsize) - pyrsize / 2, (pyrsize, 1))
467 Pangle = pup_sep * nrebin # Pyramid angle in HR pixels
468 if p_wfs.nPupils == 0:
470 # Centers is a nPupils x 2 array describing the position of the quadrant centers
471 centers = Pangle / np.sin(
472 np.pi / p_wfs.nPupils) * np.c_[np.cos(
473 (2 * np.arange(p_wfs.nPupils) + 1) * np.pi / p_wfs.nPupils),
474 np.sin((2 * np.arange(p_wfs.nPupils) + 1) *
475 np.pi / p_wfs.nPupils)]
476 # In case nPupils == 4, we put the centers back in the normal A-B-C-D ordering scheme, for misalignment processing.
477 if p_wfs.nPupils == 4:
478 centers = np.round(centers[[1, 0, 2, 3], :]).astype(np.int32)
479 # Add in the misalignments of quadrant centers, in LR px units
480 if p_wfs.pyr_misalignments is not None:
481 mis = np.asarray(p_wfs.pyr_misalignments) * nrebin
483 mis = np.zeros((p_wfs.nPupils, 2), dtype=np.float32)
484 # Pyramid as minimal intersect of 4 tilting planes
485 pyr = 2 * np.pi / pyrsize * np.min(
486 np.asarray([(c[0] + m[0]) * x + (c[1] + m[1]) * y
487 for c, m in zip(centers, mis)]), axis=0)
488 p_wfs._halfxy = np.fft.fftshift(pyr.T)
490 # Valid pixels identification
491 # Generate buffer with pupil at center
492 pup = np.zeros((pyrsize, pyrsize), dtype=np.int32)
493 pup[pyrsize // 2 - p_geom._n // 2:pyrsize // 2 + p_geom._n // 2,
494 pyrsize // 2 - p_geom._n // 2:pyrsize // 2 + p_geom._n // 2] = \
497 # Shift the mask to obtain the geometrically valid pixels
498 for qIdx in range(p_wfs.nPupils):
499 quadOnCenter = np.roll(pup, tuple((centers[qIdx] + mis[qIdx]).astype(np.int32)),
501 mskRebin = util.rebin(quadOnCenter.copy(),
502 [pyrsize // nrebin, pyrsize // nrebin])
504 stackedSubap = np.roll(mskRebin >= fracsub,
505 tuple((-centers[qIdx] / nrebin).astype(np.int32)),
507 mskRebTot = (mskRebin >= fracsub).astype(np.int32)
509 stackedSubap += np.roll(mskRebin >= fracsub,
510 tuple((-centers[qIdx] / nrebin).astype(np.int32)),
512 mskRebTot += (mskRebin >= fracsub).astype(np.int32) * (qIdx + 1)
514 validRow, validCol = [], []
515 # If n == 4, we need to tweak -again- the order to feed
516 # compass XY controllers in the appropriate order
517 if p_wfs.nPupils == 4:
518 centers = centers[[2, 1, 3, 0], :]
519 # mis = mis[[2, 1, 3, 0], :] # We don't need mis anymore - but if so keep the order of centers
520 for qIdx in range(p_wfs.nPupils):
522 np.roll(stackedSubap, tuple((centers[qIdx] / nrebin).astype(np.int32)),
524 validRow += [tmpWh[0].astype(np.int32)]
525 validCol += [tmpWh[1].astype(np.int32)]
526 nvalid = validRow[0].size
527 validRow = np.asarray(validRow).flatten()
528 validCol = np.asarray(validCol).flatten()
529 stack, index = np.unique(np.c_[validRow, validCol], axis=0, return_index=True)
531 p_wfs._nvalid = nvalid
532 p_wfs._validsubsx = validRow[np.sort(index)]
533 p_wfs._validsubsy = validCol[np.sort(index)]
534 p_wfs._hrmap = mskRebTot.astype(np.int32)
536 if (p_wfs.pyr_pos is None):
537 pixsize = (np.pi * p_wfs._qpixsize) / (3600 * 180)
538 # scale_fact = 2 * np.pi / npup * \
539 # (p_wfs.Lambda / p_tel.diam / 4.848) / pixsize * p_wfs.pyr_ampl
541 scale_fact = 2 * np.pi / pyrsize * \
542 (p_wfs.Lambda * 1e-6 / p_tel.diam) / pixsize * p_wfs.pyr_ampl
544 np.sin((np.arange(p_wfs.pyr_npts)) * 2. * np.pi / p_wfs.pyr_npts)
546 np.cos((np.arange(p_wfs.pyr_npts)) * 2. * np.pi / p_wfs.pyr_npts)
547 p_wfs.set_pyr_scale_pos(scale_fact)
548 # mod_npts = p_wfs.pyr_npts #UNUSED
551 print("Using user-defined positions [arcsec] for the pyramid modulation")
552 cx = p_wfs.pyr_pos[:, 0] / p_wfs._qpixsize
553 cy = p_wfs.pyr_pos[:, 1] / p_wfs._qpixsize
554 # mod_npts=cx.shape[0] #UNUSED
556 p_wfs._pyr_cx = cx.copy()
557 p_wfs._pyr_cy = cy.copy()
559 # telSurf = np.pi / 4. * (1 - p_tel.cobs**2.) * p_tel.diam**2.
560 # p_wfs._nphotons = p_wfs.zerop * \
561 # 10. ** (-0.4 * p_wfs.gsmag) * ittime * \
562 # p_wfs.optthroughput * telSurf
563 p_wfs._nphotons = compute_nphotons(scons.WFSType.PYRHR, ittime, p_wfs.optthroughput,
564 p_tel.diam, cobs=p_tel.cobs, zerop=p_wfs.zerop,
565 gsmag=p_wfs.gsmag, verbose=verbose)
567 # spatial filtering by the pixel extent:
568 # *2/2 intended. min should be 0.40 = sinc(0.5)^2.
569 y = np.tile(np.arange(pyrsize) - pyrsize // 2, (pyrsize, 1))
573 sincar = np.fft.fftshift(np.sinc(x) * np.sinc(y))
575 # sincar = np.roll(np.pi*x*np.pi*y,x.shape[1],axis=1)
576 p_wfs._sincar = sincar.astype(np.float32)
578 #pup = p_geom._mpupil
579 a = pyrsize // nrebin
580 b = p_geom._n // nrebin
581 pupvalid = stackedSubap[a // 2 - b // 2:a // 2 + b // 2, a // 2 - b // 2:a // 2 +
583 p_wfs._isvalid = pupvalid.T.astype(np.int32)
585 validsubsx = np.where(pupvalid)[0].astype(np.int32)
586 validsubsy = np.where(pupvalid)[1].astype(np.int32)
588 istart = np.arange(p_wfs.nxsub + 2) * p_wfs.npix
589 jstart = np.copy(istart)
591 # sorting out valid subaps
592 fluxPerSub = np.zeros((p_wfs.nxsub + 2, p_wfs.nxsub + 2), dtype=np.float32)
594 for i in range(p_wfs.nxsub + 2):
596 for j in range(p_wfs.nxsub + 2):
598 fluxPerSub[i, j] = np.sum(
599 p_geom._mpupil[indi:indi + p_wfs.npix, indj:indj + p_wfs.npix])
601 fluxPerSub = fluxPerSub / p_wfs.nxsub**2.
603 p_wfs._fluxPerSub = fluxPerSub.copy()
605 phasemap = np.zeros((p_wfs.npix * p_wfs.npix, p_wfs._nvalid), dtype=np.int32)
606 X = np.indices((p_geom._n, p_geom._n)) # we need c-like indice
607 tmp = X[1] + X[0] * p_geom._n
609 pyrtmp = np.zeros((p_geom._n, p_geom._n), dtype=np.int32)
611 for i in range(len(validsubsx)):
612 indi = istart[validsubsy[i]] # +2-1 (yorick->python
613 indj = jstart[validsubsx[i]]
614 phasemap[:, i] = tmp[indi:indi + p_wfs.npix, indj:indj + p_wfs.npix].flatten("C")
615 pyrtmp[indi:indi + p_wfs.npix, indj:indj + p_wfs.npix] = i
617 p_wfs._phasemap = phasemap
619 p_wfs._pyr_offsets = pyrtmp # pshift
def shesha.init.geom_init.init_sh_geom |
( |
conf.Param_wfs |
p_wfs, |
|
|
float |
r0, |
|
|
conf.Param_tel |
p_tel, |
|
|
conf.Param_geom |
p_geom, |
|
|
float |
ittime, |
|
|
bool |
verbose = True |
|
) |
| |
Compute the geometry of SH WFSs: valid subaps, positions of the subaps, flux per subap, etc...
:parameters: p_wfs (Param_wfs) : wfs settings
r0 (float) : atmos r0 @ 0.5 microns
p_tel (Param_tel) : telescope settings
geom (Param_geom) : geom settings
ittime (float) : 1/loop frequency [s]
verbose (bool) : (optional) display informations if True
Definition at line 639 of file geom_init.py.
642 # this is the i,j index of lower left pixel of subap in _spupil
643 istart = ((np.linspace(0, p_geom.pupdiam, p_wfs.nxsub + 1))[:-1]).astype(np.int64)
644 # Translation in _mupil useful for raytracing
646 jstart = np.copy(istart)
648 # sorting out valid subaps
649 fluxPerSub = np.zeros((p_wfs.nxsub, p_wfs.nxsub), dtype=np.float32)
651 for i in range(p_wfs.nxsub):
652 indi = istart[i] # +2-1 (yorick->python)
653 for j in range(p_wfs.nxsub):
654 indj = jstart[j] # +2-1 (yorick->python)
655 fluxPerSub[i, j] = np.sum(
656 p_geom._mpupil[indi:indi + p_wfs._pdiam, indj:indj + p_wfs._pdiam])
657 # fluxPerSub[i,j] = np.where(p_geom._mpupil[indi:indi+pdiam,indj:indj+pdiam] > 0)[0].size
659 fluxPerSub = fluxPerSub / p_wfs._pdiam**2.
661 pupvalid = (fluxPerSub >= p_wfs.fracsub) * 1
662 p_wfs._isvalid = pupvalid.astype(np.int32)
663 p_wfs._nvalid = int(np.sum(p_wfs._isvalid))
664 p_wfs._fluxPerSub = fluxPerSub.copy()
665 validx = np.where(p_wfs._isvalid.T)[1].astype(np.int32)
666 validy = np.where(p_wfs._isvalid.T)[0].astype(np.int32)
667 p_wfs._validsubsx = validx
668 p_wfs._validsubsy = validy
669 p_wfs._validpuppixx = validx * p_wfs._pdiam + 2
670 p_wfs._validpuppixy = validy * p_wfs._pdiam + 2
672 # this defines how we cut the phase into subaps
673 phasemap = np.zeros((p_wfs._pdiam * p_wfs._pdiam, p_wfs._nvalid), dtype=np.int32)
675 X = np.indices((p_geom._n, p_geom._n))
676 tmp = X[1] + X[0] * p_geom._n
679 # for i in range(p_wfs._nvalid):
681 indi = istart[p_wfs._validsubsy[i]] # +2-1 (yorick->python)
682 indj = jstart[p_wfs._validsubsx[i]]
683 phasemap[:, i] = tmp[indi:indi + p_wfs._pdiam, indj:indj +
684 p_wfs._pdiam].flatten()
685 p_wfs._phasemap = phasemap
686 p_wfs._validsubsx *= p_wfs.npix
687 p_wfs._validsubsy *= p_wfs.npix
689 # this is a phase shift of 1/2 pix in x and y
690 halfxy = np.linspace(0, 2 * np.pi, p_wfs._Nfft + 1)[0:p_wfs._pdiam] / 2.
691 halfxy = np.tile(halfxy, (p_wfs._pdiam, 1))
693 # p_wfs._halfxy = <float*>(halfxy*0.).data #dont work: half*0 is temp
696 if (p_wfs.npix % 2 == 1 and p_wfs._nrebin % 2 == 1):
697 # p_wfs._halfxy = <float*>(halfxy*0.)
698 halfxy = np.zeros((p_wfs._pdiam, p_wfs._pdiam), dtype=np.float32)
699 p_wfs._halfxy = halfxy.astype(np.float32)
701 p_wfs._halfxy = halfxy.astype(np.float32)
703 # this defines how we create a larger fov if required
704 if (p_wfs._Ntot != p_wfs._Nfft):
705 indi = int((p_wfs._Ntot - p_wfs._Nfft) / 2.) # +1 -1 (yorick>python)
706 indj = int(indi + p_wfs._Nfft)
707 X = np.indices((p_wfs._Nfft, p_wfs._Nfft)) + 1 # TODO: +1 ??
711 tmp = np.zeros((p_wfs._Ntot, p_wfs._Ntot))
712 tmp[indi:indj, indi:indj] = np.roll(x + (y - 1) * p_wfs._Nfft, p_wfs._Nfft // 2,
714 tmp[indi:indj, indi:indj] = np.roll(tmp[indi:indj, indi:indj], p_wfs._Nfft // 2,
717 tmp = np.roll(tmp, p_wfs._Ntot // 2, axis=0)
718 tmp = np.roll(tmp, p_wfs._Ntot // 2, axis=1)
720 tmp = np.where(tmp.flatten())[0]
722 p_wfs._hrmap = np.copy(tmp.astype(np.int32))
723 p_wfs._hrmap = np.reshape(p_wfs._hrmap, (p_wfs._hrmap.shape[0], 1))
724 # must be set even if unused
727 tmp = np.zeros((2, 2))
728 p_wfs._hrmap = np.copy(tmp.astype(np.int32))
729 # must be set even if unused
731 if (p_wfs._nrebin * p_wfs.npix % 2 != p_wfs._Ntot % 2):
733 indi = int((p_wfs._Ntot - p_wfs._nrebin * p_wfs.npix) / 2.) + 1
735 indi = int((p_wfs._Ntot - p_wfs._nrebin * p_wfs.npix) / 2.) + 0
736 indj = int(indi + p_wfs._nrebin * p_wfs.npix)
738 X = np.indices((p_wfs._nrebin * p_wfs.npix, p_wfs._nrebin * p_wfs.npix))
739 x = (X[1] / p_wfs._nrebin).astype(np.int64)
740 y = (X[0] / p_wfs._nrebin).astype(np.int64)
743 binindices = np.zeros((p_wfs._Ntot, p_wfs._Ntot))
744 binindices[indi:indj, indi:indj] = x + y * p_wfs.npix + 1
746 binmap = np.zeros((p_wfs._nrebin * p_wfs._nrebin, p_wfs.npix * p_wfs.npix))
748 X = np.indices((p_wfs._Ntot, p_wfs._Ntot))
749 tmp = X[1] + X[0] * p_wfs._Ntot
751 if (p_wfs.gsalt <= 0):
752 binindices = np.roll(binindices, binindices.shape[0] // 2, axis=0)
753 binindices = np.roll(binindices, binindices.shape[1] // 2, axis=1)
755 for i in range(p_wfs.npix * p_wfs.npix):
756 binmap[:, i] = tmp[np.where(binindices == i + 1)]
757 # binmap=np.reshape(binmap.flatten("F"),(binmap.shape[0],binmap.shape[1]))
758 p_wfs._binmap = np.copy(binmap.astype(np.int32))
760 dr0 = p_tel.diam / r0 * \
761 (0.5 / p_wfs.Lambda) ** 1.2 / \
762 np.cos(p_geom.zenithangle * CONST.DEG2RAD) ** 0.6
763 fwhmseeing = p_wfs.Lambda / \
764 (p_tel.diam / np.sqrt(p_wfs.nxsub ** 2. + (dr0 / 1.5) ** 2.)) / 4.848
765 kernelfwhm = np.sqrt(fwhmseeing**2. + p_wfs.kernel**2.)
767 tmp = util.makegaussian(p_wfs._Ntot, kernelfwhm / p_wfs._qpixsize, p_wfs._Ntot // 2,
768 p_wfs._Ntot // 2).astype(np.float32)
770 tmp = np.roll(tmp, tmp.shape[0] // 2, axis=0)
771 tmp = np.roll(tmp, tmp.shape[1] // 2, axis=1)
773 tmp[0, 0] = 1. # this insures that even with fwhm=0, the kernel is a dirac
774 tmp = tmp / np.sum(tmp)
775 tmp = np.fft.fft2(tmp).astype(np.complex64) / (p_wfs._Ntot * p_wfs._Ntot)
776 p_wfs._ftkernel = np.copy(tmp).astype(np.complex64)
778 # dealing with photometry
779 # telSurf = np.pi / 4. * (1 - p_tel.cobs**2.) * p_tel.diam**2.
781 # from the guide star
782 if (p_wfs.gsalt == 0):
783 if (p_wfs.zerop == 0):
785 # p_wfs._nphotons = p_wfs.zerop * 10 ** (-0.4 * p_wfs.gsmag) * \
786 # p_wfs.optthroughput * \
787 # (p_tel.diam / p_wfs.nxsub) ** 2. * ittime
788 # include throughput to WFS
789 # for unobstructed subaperture
791 p_wfs._nphotons = compute_nphotons(scons.WFSType.SH, ittime, p_wfs.optthroughput,
792 p_tel.diam, nxsub=p_wfs.nxsub,
793 zerop=p_wfs.zerop, gsmag=p_wfs.gsmag,
796 else: # we are dealing with a LGS
797 # p_wfs._nphotons = p_wfs.lgsreturnperwatt * \
798 # p_wfs.laserpower * \
799 # p_wfs.optthroughput * \
800 # (p_tel.diam / p_wfs.nxsub) ** 2. * 1e4 * \
803 # ... for given power
804 # include throughput to WFS
805 # for unobstructed subaperture
807 p_wfs._nphotons = compute_nphotons(scons.WFSType.SH, ittime, p_wfs.optthroughput,
808 p_tel.diam, nxsub=p_wfs.nxsub,
809 lgsreturnperwatt=p_wfs.lgsreturnperwatt,
810 laserpower=p_wfs.laserpower, verbose=verbose)
def shesha.init.geom_init.init_wfs_geom |
( |
conf.Param_wfs |
p_wfs, |
|
|
float |
r0, |
|
|
conf.Param_tel |
p_tel, |
|
|
conf.Param_geom |
p_geom, |
|
|
float |
ittime, |
|
|
|
verbose = 1 |
|
) |
| |
Compute the geometry of WFSs: valid subaps, positions of the subaps, flux per subap, etc...
:parameters: p_wfs (Param_wfs) : wfs settings
r0 (float) : atmos r0 @ 0.5 microns
p_tel (Param_tel) : telescope settings
geom (Param_geom) : geom settings
ittime (float) : 1/loop frequency [s]
verbose (int) : (optional) display informations if 0
Definition at line 129 of file geom_init.py.
133 if p_wfs.type == scons.WFSType.SH or p_wfs.type == scons.WFSType.PYRHR or p_wfs.type == scons.WFSType.PYRLR:
134 pdiam = p_geom.pupdiam // p_wfs.nxsub
135 if (p_geom.pupdiam % p_wfs.nxsub > 0):
143 init_wfs_size(p_wfs, r0, p_tel, verbose)
145 if not p_geom.is_init:
146 # this is the wfs with largest # of subaps
147 # the overall geometry is deduced from it
148 #if not p_geom.pupdiam:
149 if p_geom.pupdiam != 0 and p_geom.pupdiam != p_wfs._pdiam * p_wfs.nxsub:
150 print("WARNING: pupdiam set value not correct")
151 p_geom.pupdiam = p_wfs._pdiam * p_wfs.nxsub
152 print("pupdiam used: ", p_geom.pupdiam)
153 if p_wfs.type == scons.WFSType.PYRHR or p_wfs.type == scons.WFSType.PYRLR:
154 geom_init(p_geom, p_tel, padding=p_wfs._nrebin)
155 elif (p_wfs.type == scons.WFSType.SH):
156 geom_init(p_geom, p_tel)
158 raise RuntimeError("This WFS can not be used")
160 if (p_wfs.type == scons.WFSType.PYRHR or p_wfs.type == scons.WFSType.PYRLR):
161 init_pyrhr_geom(p_wfs, r0, p_tel, p_geom, ittime, verbose=1)
162 elif (p_wfs.type == scons.WFSType.SH):
163 init_sh_geom(p_wfs, r0, p_tel, p_geom, ittime, verbose=1)
165 raise RuntimeError("This WFS can not be used")
def shesha.init.geom_init.init_wfs_size |
( |
conf.Param_wfs |
p_wfs, |
|
|
float |
r0, |
|
|
conf.Param_tel |
p_tel, |
|
|
|
verbose = 1 |
|
) |
| |
Compute all the parameters usefull for further WFS image computation (array sizes)
:parameters: p_wfs (Param_wfs) : wfs settings
r0 (float) : atmos r0 @ 0.5 microns
p_tel (Param_tel) : telescope settings
verbose (int) : (optional) display informations if 0
Scheme to determine arrays sizes sh : k = 6 p = k * d/r0 n = int(2*d*v/lambda/CONST.RAD2ARCSEC)+1 N = fft_goodsize(k*n/v*lambda/r0*CONST.RAD2ARCSEC) u = k * lambda / r0 * CONST.RAD2ARCSEC / N n = v/u - int(v/u) > 0.5 ? int(v/u)+1 : int(v/u) v = n * u Nt = v * Npix
pyr : Fs = field stop radius in arcsec N size of big array for FFT in pixels P pupil diameter in pixels D diameter of telescope in m Nssp number of pyr measurement points in the pupil
Rf = Fs . N . D / lambda / P ideally we choose : Fs = lambda / D . Nssp / 2
if we want good sampling of r0 (avoid aliasing of speckles) P > D / r0 . m with m = 2 or 3
to get reasonable space between pupil images : N > P.(2 + 3S) with S close to 1 N must be a power of 2
to ease computation lets assume that Nssp is a divider of P scaling factor between high res pupil images in pyramid model and actual pupil size on camera images would be P / Nssp
Definition at line 213 of file geom_init.py.
215 r0 = r0 * (p_wfs.Lambda * 2)**(6. / 5)
219 print("r0 for WFS :", "%3.2f" % r0, " m")
220 # seeing = CONST.RAD2ARCSEC * (p_wfs.lambda * 1.e-6) / r0
222 print("seeing for WFS : ",
223 "%3.2f" % (CONST.RAD2ARCSEC * (p_wfs.Lambda * 1.e-6) / r0), "\"")
225 if (p_wfs._pdiam <= 0):
226 # this case is usualy for the wfs with max # of subaps
227 # we look for the best compromise between pixsize and fov
228 if (p_wfs.type == scons.WFSType.SH):
230 subapdiam = p_tel.diam / float(p_wfs.nxsub) # diam of subap
232 pdiam = int(k * subapdiam / r0) # number of phase points per subap
236 # Must be even to keep ssp and actuators grids aligned in the pupil
237 if ((pdiam * p_wfs.nxsub) % 2):
240 nrebin = int(2 * subapdiam * p_wfs.pixsize /
241 (p_wfs.Lambda * 1.e-6) / CONST.RAD2ARCSEC) + 1
242 nrebin = max(2, nrebin)
243 # first atempt on a rebin factor
245 # since we clipped pdiam we have to be carreful in nfft computation
246 Nfft = util.fft_goodsize(
247 int(pdiam / subapdiam * nrebin / p_wfs.pixsize * CONST.RAD2ARCSEC *
248 (p_wfs.Lambda * 1.e-6)))
250 elif (p_wfs.type == scons.WFSType.PYRHR or p_wfs.type == scons.WFSType.PYRLR):
251 # while (pdiam % p_wfs.npix != 0) pdiam+=1;
252 k = 3 if p_wfs.type == scons.WFSType.PYRHR else 2
253 pdiam = int(p_tel.diam / r0 * k)
254 while (pdiam % p_wfs.nxsub != 0):
255 pdiam += 1 # we choose to have a multiple of p_wfs.nxsub
256 pdiam = pdiam // p_wfs.nxsub
263 # this case is for a wfs with fixed # of phase points
264 Nfft = util.fft_goodsize(2 * pdiam)
265 # size of the support in fourier domain
267 # qpixsize = pdiam * \
268 # (p_wfs.Lambda * 1.e-6) / subapdiam * CONST.RAD2ARCSEC / Nfft
269 # # quantum pixel size
271 if (p_wfs.type == scons.WFSType.SH):
272 subapdiam = p_tel.diam / float(p_wfs.nxsub) # diam of subap
274 # size of the support in fourier domain
276 # qpixsize = k * (p_wfs.Lambda*1.e-6) / r0 * CONST.RAD2ARCSEC / Nfft
277 qpixsize = (pdiam * (p_wfs.Lambda * 1.e-6) / subapdiam * CONST.RAD2ARCSEC) / Nfft
279 # # actual rebin factor
280 if (p_wfs.pixsize / qpixsize - int(p_wfs.pixsize / qpixsize) > 0.5):
281 nrebin = int(p_wfs.pixsize / qpixsize) + 1
283 nrebin = int(p_wfs.pixsize / qpixsize)
286 pixsize = nrebin * qpixsize
288 if (pixsize * p_wfs.npix > qpixsize * Nfft):
289 Ntot = util.fft_goodsize(int(pixsize * p_wfs.npix / qpixsize) + 1)
293 if (Ntot % 2 != Nfft % 2):
296 elif (p_wfs.type == scons.WFSType.PYRHR or p_wfs.type == scons.WFSType.PYRLR):
297 pdiam = pdiam * p_wfs.nxsub
299 # fft_goodsize( m * pdiam)
300 Nfft = int(2**np.ceil(np.log2(m * pdiam)))
302 nrebin = pdiam // p_wfs.nxsub
303 while (Nfft % nrebin != 0):
304 nrebin += 1 # we choose to have a divisor of Nfft
305 pdiam = nrebin * p_wfs.nxsub
306 Nfft = int(2**np.ceil(np.log2(m * pdiam)))
309 (p_wfs.Lambda * 1.e-6) / p_tel.diam * CONST.RAD2ARCSEC) / Nfft
311 Ntot = Nfft // pdiam * p_wfs.nxsub
312 pixsize = qpixsize * nrebin
313 pdiam = pdiam // p_wfs.nxsub
316 p_wfs.pixsize = pixsize
317 p_wfs._qpixsize = qpixsize
320 p_wfs._nrebin = nrebin
321 p_wfs._subapd = p_tel.diam / p_wfs.nxsub
324 if (p_wfs.type == scons.WFSType.SH):
325 print("quantum pixsize : ", "%5.4f" % qpixsize, "\"")
326 print("simulated FoV : ", "%3.2f" % (Ntot * qpixsize), "\" x ",
327 "%3.2f" % (Ntot * qpixsize), "\"")
328 print("actual pixsize : ", "%5.4f" % pixsize)
329 print("actual FoV : ", "%3.2f" % (pixsize * p_wfs.npix), "\" x ",
330 "%3.2f" % (pixsize * p_wfs.npix), "\"")
331 print("number of phase points : ", p_wfs._pdiam)
332 print("size of fft support : ", Nfft)
333 print("size of HR spot support : ", Ntot)
335 elif (p_wfs.type == scons.WFSType.PYRHR or p_wfs.type == scons.WFSType.PYRLR):
336 print("quantum pixsize in pyr image : ", "%5.4f" % qpixsize, "\"")
337 print("simulated FoV : ", "%3.2f" % (Nfft * qpixsize), "\" x ",
338 "%3.2f" % (Nfft * qpixsize), "\"")
339 print("number of phase points : ", p_wfs._pdiam * p_wfs.nxsub)
340 print("size of fft support : ", Nfft)