LSST Applications g0265f82a02+c6dfa2ddaf,g1162b98a3f+b2075782a9,g2079a07aa2+1b2e822518,g2bbee38e9b+c6dfa2ddaf,g337abbeb29+c6dfa2ddaf,g3ddfee87b4+a60788ef87,g50ff169b8f+2eb0e556e8,g52b1c1532d+90ebb246c7,g555ede804d+a60788ef87,g591dd9f2cf+ba8caea58f,g5ec818987f+864ee9cddb,g858d7b2824+9ee1ab4172,g876c692160+a40945ebb7,g8a8a8dda67+90ebb246c7,g8cdfe0ae6a+4fd9e222a8,g99cad8db69+5e309b7bc6,g9ddcbc5298+a1346535a5,ga1e77700b3+df8f93165b,ga8c6da7877+aa12a14d27,gae46bcf261+c6dfa2ddaf,gb0e22166c9+8634eb87fb,gb3f2274832+d0da15e3be,gba4ed39666+1ac82b564f,gbb8dafda3b+5dfd9c994b,gbeb006f7da+97157f9740,gc28159a63d+c6dfa2ddaf,gc86a011abf+9ee1ab4172,gcf0d15dbbd+a60788ef87,gdaeeff99f8+1cafcb7cd4,gdc0c513512+9ee1ab4172,ge79ae78c31+c6dfa2ddaf,geb67518f79+ba1859f325,geb961e4c1e+f9439d1e6f,gee10cc3b42+90ebb246c7,gf1cff7945b+9ee1ab4172,w.2024.12
LSST Data Management Base Package
Loading...
Searching...
No Matches
Public Member Functions | Public Attributes | Static Public Attributes | Static Protected Attributes | List of all members
lsst.ip.isr.isrMock.IsrMock Class Reference
Inheritance diagram for lsst.ip.isr.isrMock.IsrMock:
lsst.ip.isr.isrMock.BfKernelMock lsst.ip.isr.isrMock.CrosstalkCoeffMock lsst.ip.isr.isrMock.DefectMock lsst.ip.isr.isrMock.MasterMock lsst.ip.isr.isrMock.RawMock lsst.ip.isr.isrMock.TransmissionMock lsst.ip.isr.isrMock.BiasMock lsst.ip.isr.isrMock.DarkMock lsst.ip.isr.isrMock.FlatMock lsst.ip.isr.isrMock.FringeMock lsst.ip.isr.isrMock.CalibratedRawMock lsst.ip.isr.isrMock.RawDictMock lsst.ip.isr.isrMock.TrimmedRawMock lsst.ip.isr.isrMock.UntrimmedFringeMock

Public Member Functions

 __init__ (self, **kwargs)
 
 run (self)
 
 makeData (self)
 
 makeBfKernel (self)
 
 makeDefectList (self)
 
 makeCrosstalkCoeff (self)
 
 makeTransmissionCurve (self)
 
 makeImage (self)
 
 getCamera (self)
 
 getExposure (self)
 
 getWcs (self)
 
 localCoordToExpCoord (self, ampData, x, y)
 
 amplifierAddNoise (self, ampData, mean, sigma)
 
 amplifierAddYGradient (self, ampData, start, end)
 
 amplifierAddSource (self, ampData, scale, x0, y0)
 
 amplifierAddCT (self, ampDataSource, ampDataTarget, scale)
 
 amplifierAddFringe (self, amp, ampData, scale, x0=100, y0=0)
 
 amplifierMultiplyFlat (self, amp, ampData, fracDrop, u0=100.0, v0=100.0)
 

Public Attributes

 rng
 
 crosstalkCoeffs
 
 bfKernel
 

Static Public Attributes

 ConfigClass = IsrMockConfig
 

Static Protected Attributes

str _DefaultName = "isrMock"
 

Detailed Description

Class to generate consistent mock images for ISR testing.

ISR testing currently relies on one-off fake images that do not
accurately mimic the full set of detector effects. This class
uses the test camera/detector/amplifier structure defined in
`lsst.afw.cameraGeom.testUtils` to avoid making the test data
dependent on any of the actual obs package formats.

Definition at line 253 of file isrMock.py.

Constructor & Destructor Documentation

◆ __init__()

lsst.ip.isr.isrMock.IsrMock.__init__ ( self,
** kwargs )

Reimplemented in lsst.ip.isr.isrMock.RawMock, lsst.ip.isr.isrMock.TrimmedRawMock, lsst.ip.isr.isrMock.CalibratedRawMock, lsst.ip.isr.isrMock.RawDictMock, lsst.ip.isr.isrMock.MasterMock, lsst.ip.isr.isrMock.BiasMock, lsst.ip.isr.isrMock.DarkMock, lsst.ip.isr.isrMock.FlatMock, lsst.ip.isr.isrMock.FringeMock, lsst.ip.isr.isrMock.UntrimmedFringeMock, lsst.ip.isr.isrMock.BfKernelMock, lsst.ip.isr.isrMock.DefectMock, lsst.ip.isr.isrMock.CrosstalkCoeffMock, and lsst.ip.isr.isrMock.TransmissionMock.

Definition at line 265 of file isrMock.py.

265 def __init__(self, **kwargs):
266 super().__init__(**kwargs)
267 self.rng = np.random.RandomState(self.config.rngSeed)
268 self.crosstalkCoeffs = np.array([[0.0, 0.0, 0.0, 0.0, 0.0, -1e-3, 0.0, 0.0],
269 [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
270 [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
271 [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
272 [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
273 [1e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
274 [1e-2, 0.0, 0.0, 2.2e-2, 0.0, 0.0, 0.0, 0.0],
275 [1e-2, 5e-3, 5e-4, 3e-3, 4e-2, 5e-3, 5e-3, 0.0]])
276
277 self.bfKernel = np.array([[1., 4., 7., 4., 1.],
278 [4., 16., 26., 16., 4.],
279 [7., 26., 41., 26., 7.],
280 [4., 16., 26., 16., 4.],
281 [1., 4., 7., 4., 1.]]) / 273.0
282

Member Function Documentation

◆ amplifierAddCT()

lsst.ip.isr.isrMock.IsrMock.amplifierAddCT ( self,
ampDataSource,
ampDataTarget,
scale )
Add a scaled copy of an amplifier to another, simulating crosstalk.

 This method operates in the amplifier coordinate frame.

Parameters
----------
ampDataSource : `lsst.afw.image.ImageF`
    Amplifier image to add scaled copy from.
ampDataTarget : `lsst.afw.image.ImageF`
    Amplifier image to add scaled copy to.
scale : `float`
    Flux scale of the copy to add to the target.

Notes
-----
This simulates simple crosstalk between amplifiers.

Definition at line 681 of file isrMock.py.

681 def amplifierAddCT(self, ampDataSource, ampDataTarget, scale):
682 """Add a scaled copy of an amplifier to another, simulating crosstalk.
683
684 This method operates in the amplifier coordinate frame.
685
686 Parameters
687 ----------
688 ampDataSource : `lsst.afw.image.ImageF`
689 Amplifier image to add scaled copy from.
690 ampDataTarget : `lsst.afw.image.ImageF`
691 Amplifier image to add scaled copy to.
692 scale : `float`
693 Flux scale of the copy to add to the target.
694
695 Notes
696 -----
697 This simulates simple crosstalk between amplifiers.
698 """
699 ampDataTarget.array[:] = (ampDataTarget.array[:]
700 + scale * ampDataSource.array[:])
701

◆ amplifierAddFringe()

lsst.ip.isr.isrMock.IsrMock.amplifierAddFringe ( self,
amp,
ampData,
scale,
x0 = 100,
y0 = 0 )
Add a fringe-like ripple pattern to an amplifier's image data.

Parameters
----------
amp : `~lsst.afw.ampInfo.AmpInfoRecord`
    Amplifier to operate on. Needed for amp<->exp coordinate
    transforms.
ampData : `lsst.afw.image.ImageF`
    Amplifier image to operate on.
scale : `numpy.array` or `float`
    Peak intensity scaling for the ripple.
x0 : `numpy.array` or `float`, optional
    Fringe center
y0 : `numpy.array` or `float`, optional
    Fringe center

Notes
-----
This uses an offset sinc function to generate a ripple
pattern. True fringes have much finer structure, but this
pattern should be visually identifiable. The (x, y)
coordinates are in the frame of the amplifier, and (u, v) in
the frame of the full trimmed image.

Definition at line 703 of file isrMock.py.

703 def amplifierAddFringe(self, amp, ampData, scale, x0=100, y0=0):
704 """Add a fringe-like ripple pattern to an amplifier's image data.
705
706 Parameters
707 ----------
708 amp : `~lsst.afw.ampInfo.AmpInfoRecord`
709 Amplifier to operate on. Needed for amp<->exp coordinate
710 transforms.
711 ampData : `lsst.afw.image.ImageF`
712 Amplifier image to operate on.
713 scale : `numpy.array` or `float`
714 Peak intensity scaling for the ripple.
715 x0 : `numpy.array` or `float`, optional
716 Fringe center
717 y0 : `numpy.array` or `float`, optional
718 Fringe center
719
720 Notes
721 -----
722 This uses an offset sinc function to generate a ripple
723 pattern. True fringes have much finer structure, but this
724 pattern should be visually identifiable. The (x, y)
725 coordinates are in the frame of the amplifier, and (u, v) in
726 the frame of the full trimmed image.
727 """
728 for x in range(0, ampData.getDimensions().getX()):
729 for y in range(0, ampData.getDimensions().getY()):
730 (u, v) = self.localCoordToExpCoord(amp, x, y)
731 ampData.getArray()[y][x] = np.sum((ampData.getArray()[y][x]
732 + scale * np.sinc(((u - x0) / 50)**2
733 + ((v - y0) / 50)**2)))
734

◆ amplifierAddNoise()

lsst.ip.isr.isrMock.IsrMock.amplifierAddNoise ( self,
ampData,
mean,
sigma )
Add Gaussian noise to an amplifier's image data.

 This method operates in the amplifier coordinate frame.

Parameters
----------
ampData : `lsst.afw.image.ImageF`
    Amplifier image to operate on.
mean : `float`
    Mean value of the Gaussian noise.
sigma : `float`
    Sigma of the Gaussian noise.

Definition at line 623 of file isrMock.py.

623 def amplifierAddNoise(self, ampData, mean, sigma):
624 """Add Gaussian noise to an amplifier's image data.
625
626 This method operates in the amplifier coordinate frame.
627
628 Parameters
629 ----------
630 ampData : `lsst.afw.image.ImageF`
631 Amplifier image to operate on.
632 mean : `float`
633 Mean value of the Gaussian noise.
634 sigma : `float`
635 Sigma of the Gaussian noise.
636 """
637 ampArr = ampData.array
638 ampArr[:] = ampArr[:] + self.rng.normal(mean, sigma,
639 size=ampData.getDimensions()).transpose()
640

◆ amplifierAddSource()

lsst.ip.isr.isrMock.IsrMock.amplifierAddSource ( self,
ampData,
scale,
x0,
y0 )
Add a single Gaussian source to an amplifier.

 This method operates in the amplifier coordinate frame.

Parameters
----------
ampData : `lsst.afw.image.ImageF`
    Amplifier image to operate on.
scale : `float`
    Peak flux of the source to add.
x0 : `float`
    X-coordinate of the source peak.
y0 : `float`
    Y-coordinate of the source peak.

Definition at line 660 of file isrMock.py.

660 def amplifierAddSource(self, ampData, scale, x0, y0):
661 """Add a single Gaussian source to an amplifier.
662
663 This method operates in the amplifier coordinate frame.
664
665 Parameters
666 ----------
667 ampData : `lsst.afw.image.ImageF`
668 Amplifier image to operate on.
669 scale : `float`
670 Peak flux of the source to add.
671 x0 : `float`
672 X-coordinate of the source peak.
673 y0 : `float`
674 Y-coordinate of the source peak.
675 """
676 for x in range(0, ampData.getDimensions().getX()):
677 for y in range(0, ampData.getDimensions().getY()):
678 ampData.array[y][x] = (ampData.array[y][x]
679 + scale * np.exp(-0.5 * ((x - x0)**2 + (y - y0)**2) / 3.0**2))
680

◆ amplifierAddYGradient()

lsst.ip.isr.isrMock.IsrMock.amplifierAddYGradient ( self,
ampData,
start,
end )
Add a y-axis linear gradient to an amplifier's image data.

 This method operates in the amplifier coordinate frame.

Parameters
----------
ampData : `lsst.afw.image.ImageF`
    Amplifier image to operate on.
start : `float`
    Start value of the gradient (at y=0).
end : `float`
    End value of the gradient (at y=ymax).

Definition at line 641 of file isrMock.py.

641 def amplifierAddYGradient(self, ampData, start, end):
642 """Add a y-axis linear gradient to an amplifier's image data.
643
644 This method operates in the amplifier coordinate frame.
645
646 Parameters
647 ----------
648 ampData : `lsst.afw.image.ImageF`
649 Amplifier image to operate on.
650 start : `float`
651 Start value of the gradient (at y=0).
652 end : `float`
653 End value of the gradient (at y=ymax).
654 """
655 nPixY = ampData.getDimensions().getY()
656 ampArr = ampData.array
657 ampArr[:] = ampArr[:] + (np.interp(range(nPixY), (0, nPixY - 1), (start, end)).reshape(nPixY, 1)
658 + np.zeros(ampData.getDimensions()).transpose())
659

◆ amplifierMultiplyFlat()

lsst.ip.isr.isrMock.IsrMock.amplifierMultiplyFlat ( self,
amp,
ampData,
fracDrop,
u0 = 100.0,
v0 = 100.0 )
Multiply an amplifier's image data by a flat-like pattern.

Parameters
----------
amp : `lsst.afw.ampInfo.AmpInfoRecord`
    Amplifier to operate on. Needed for amp<->exp coordinate
    transforms.
ampData : `lsst.afw.image.ImageF`
    Amplifier image to operate on.
fracDrop : `float`
    Fractional drop from center to edge of detector along x-axis.
u0 : `float`
    Peak location in detector coordinates.
v0 : `float`
    Peak location in detector coordinates.

Notes
-----
This uses a 2-d Gaussian to simulate an illumination pattern
that falls off towards the edge of the detector. The (x, y)
coordinates are in the frame of the amplifier, and (u, v) in
the frame of the full trimmed image.

Definition at line 735 of file isrMock.py.

735 def amplifierMultiplyFlat(self, amp, ampData, fracDrop, u0=100.0, v0=100.0):
736 """Multiply an amplifier's image data by a flat-like pattern.
737
738 Parameters
739 ----------
740 amp : `lsst.afw.ampInfo.AmpInfoRecord`
741 Amplifier to operate on. Needed for amp<->exp coordinate
742 transforms.
743 ampData : `lsst.afw.image.ImageF`
744 Amplifier image to operate on.
745 fracDrop : `float`
746 Fractional drop from center to edge of detector along x-axis.
747 u0 : `float`
748 Peak location in detector coordinates.
749 v0 : `float`
750 Peak location in detector coordinates.
751
752 Notes
753 -----
754 This uses a 2-d Gaussian to simulate an illumination pattern
755 that falls off towards the edge of the detector. The (x, y)
756 coordinates are in the frame of the amplifier, and (u, v) in
757 the frame of the full trimmed image.
758 """
759 if fracDrop >= 1.0:
760 raise RuntimeError("Flat fractional drop cannot be greater than 1.0")
761
762 sigma = u0 / np.sqrt(-2.0 * np.log(fracDrop))
763
764 for x in range(0, ampData.getDimensions().getX()):
765 for y in range(0, ampData.getDimensions().getY()):
766 (u, v) = self.localCoordToExpCoord(amp, x, y)
767 f = np.exp(-0.5 * ((u - u0)**2 + (v - v0)**2) / sigma**2)
768 ampData.array[y][x] = (ampData.array[y][x] * f)
769
770

◆ getCamera()

lsst.ip.isr.isrMock.IsrMock.getCamera ( self)
Construct a test camera object.

Returns
-------
camera : `lsst.afw.cameraGeom.camera`
    Test camera.

Definition at line 501 of file isrMock.py.

501 def getCamera(self):
502 """Construct a test camera object.
503
504 Returns
505 -------
506 camera : `lsst.afw.cameraGeom.camera`
507 Test camera.
508 """
509 cameraWrapper = afwTestUtils.CameraWrapper(
510 plateScale=self.config.plateScale,
511 radialDistortion=self.config.radialDistortion,
512 isLsstLike=self.config.isLsstLike,
513 )
514 camera = cameraWrapper.camera
515 return camera
516

◆ getExposure()

lsst.ip.isr.isrMock.IsrMock.getExposure ( self)
Construct a test exposure.

The test exposure has a simple WCS set, as well as a list of
unlikely header keywords that can be removed during ISR
processing to exercise that code.

Returns
-------
exposure : `lsst.afw.exposure.Exposure`
    Construct exposure containing masked image of the
    appropriate size.

Definition at line 517 of file isrMock.py.

517 def getExposure(self):
518 """Construct a test exposure.
519
520 The test exposure has a simple WCS set, as well as a list of
521 unlikely header keywords that can be removed during ISR
522 processing to exercise that code.
523
524 Returns
525 -------
526 exposure : `lsst.afw.exposure.Exposure`
527 Construct exposure containing masked image of the
528 appropriate size.
529 """
530 camera = self.getCamera()
531 detector = camera[self.config.detectorIndex]
532 image = afwUtils.makeImageFromCcd(detector,
533 isTrimmed=self.config.isTrimmed,
534 showAmpGain=False,
535 rcMarkSize=0,
536 binSize=1,
537 imageFactory=afwImage.ImageF)
538
539 var = afwImage.ImageF(image.getDimensions())
540 mask = afwImage.Mask(image.getDimensions())
541 image.assign(0.0)
542
543 maskedImage = afwImage.makeMaskedImage(image, mask, var)
544 exposure = afwImage.makeExposure(maskedImage)
545 exposure.setDetector(detector)
546 exposure.setWcs(self.getWcs())
547
548 visitInfo = afwImage.VisitInfo(exposureTime=self.config.expTime, darkTime=self.config.darkTime)
549 exposure.getInfo().setVisitInfo(visitInfo)
550
551 metadata = exposure.getMetadata()
552 metadata.add("SHEEP", 7.3, "number of sheep on farm")
553 metadata.add("MONKEYS", 155, "monkeys per tree")
554 metadata.add("VAMPIRES", 4, "How scary are vampires.")
555
556 ccd = exposure.getDetector()
557 newCcd = ccd.rebuild()
558 newCcd.clear()
559 for amp in ccd:
560 newAmp = amp.rebuild()
561 newAmp.setLinearityCoeffs((0., 1., 0., 0.))
562 newAmp.setLinearityType("Polynomial")
563 newAmp.setGain(self.config.gain)
564 newAmp.setSuspectLevel(25000.0)
565 newAmp.setSaturation(32000.0)
566 newCcd.append(newAmp)
567 exposure.setDetector(newCcd.finish())
568
569 exposure.image.array[:] = np.zeros(exposure.getImage().getDimensions()).transpose()
570 exposure.mask.array[:] = np.zeros(exposure.getMask().getDimensions()).transpose()
571 exposure.variance.array[:] = np.zeros(exposure.getVariance().getDimensions()).transpose()
572
573 return exposure
574
Represent a 2-dimensional array of bitmask pixels.
Definition Mask.h:77
Information about a single exposure of an imaging camera.
Definition VisitInfo.h:68
MaskedImage< ImagePixelT, MaskPixelT, VariancePixelT > * makeMaskedImage(typename std::shared_ptr< Image< ImagePixelT > > image, typename std::shared_ptr< Mask< MaskPixelT > > mask=Mask< MaskPixelT >(), typename std::shared_ptr< Image< VariancePixelT > > variance=Image< VariancePixelT >())
A function to return a MaskedImage of the correct type (cf.
std::shared_ptr< Exposure< ImagePixelT, MaskPixelT, VariancePixelT > > makeExposure(MaskedImage< ImagePixelT, MaskPixelT, VariancePixelT > &mimage, std::shared_ptr< geom::SkyWcs const > wcs=std::shared_ptr< geom::SkyWcs const >())
A function to return an Exposure of the correct type (cf.
Definition Exposure.h:484

◆ getWcs()

lsst.ip.isr.isrMock.IsrMock.getWcs ( self)
Construct a dummy WCS object.

Taken from the deprecated ip_isr/examples/exampleUtils.py.

This is not guaranteed, given the distortion and pixel scale
listed in the afwTestUtils camera definition.

Returns
-------
wcs : `lsst.afw.geom.SkyWcs`
    Test WCS transform.

Definition at line 575 of file isrMock.py.

575 def getWcs(self):
576 """Construct a dummy WCS object.
577
578 Taken from the deprecated ip_isr/examples/exampleUtils.py.
579
580 This is not guaranteed, given the distortion and pixel scale
581 listed in the afwTestUtils camera definition.
582
583 Returns
584 -------
585 wcs : `lsst.afw.geom.SkyWcs`
586 Test WCS transform.
587 """
588 return afwGeom.makeSkyWcs(crpix=lsst.geom.Point2D(0.0, 100.0),
589 crval=lsst.geom.SpherePoint(45.0, 25.0, lsst.geom.degrees),
590 cdMatrix=afwGeom.makeCdMatrix(scale=1.0*lsst.geom.degrees))
591
Point in an unspecified spherical coordinate system.
Definition SpherePoint.h:57
std::shared_ptr< SkyWcs > makeSkyWcs(daf::base::PropertySet &metadata, bool strip=false)
Construct a SkyWcs from FITS keywords.
Definition SkyWcs.cc:521
Eigen::Matrix2d makeCdMatrix(lsst::geom::Angle const &scale, lsst::geom::Angle const &orientation=0 *lsst::geom::degrees, bool flipX=false)
Make a WCS CD matrix.
Definition SkyWcs.cc:133

◆ localCoordToExpCoord()

lsst.ip.isr.isrMock.IsrMock.localCoordToExpCoord ( self,
ampData,
x,
y )
Convert between a local amplifier coordinate and the full
exposure coordinate.

Parameters
----------
ampData : `lsst.afw.image.ImageF`
    Amplifier image to use for conversions.
x : `int`
    X-coordinate of the point to transform.
y : `int`
    Y-coordinate of the point to transform.

Returns
-------
u : `int`
    Transformed x-coordinate.
v : `int`
    Transformed y-coordinate.

Notes
-----
The output is transposed intentionally here, to match the
internal transpose between numpy and afw.image coordinates.

Definition at line 592 of file isrMock.py.

592 def localCoordToExpCoord(self, ampData, x, y):
593 """Convert between a local amplifier coordinate and the full
594 exposure coordinate.
595
596 Parameters
597 ----------
598 ampData : `lsst.afw.image.ImageF`
599 Amplifier image to use for conversions.
600 x : `int`
601 X-coordinate of the point to transform.
602 y : `int`
603 Y-coordinate of the point to transform.
604
605 Returns
606 -------
607 u : `int`
608 Transformed x-coordinate.
609 v : `int`
610 Transformed y-coordinate.
611
612 Notes
613 -----
614 The output is transposed intentionally here, to match the
615 internal transpose between numpy and afw.image coordinates.
616 """
617 u = x + ampData.getBBox().getBeginX()
618 v = y + ampData.getBBox().getBeginY()
619
620 return (v, u)
621

◆ makeBfKernel()

lsst.ip.isr.isrMock.IsrMock.makeBfKernel ( self)
Generate a simple Gaussian brighter-fatter kernel.

Returns
-------
kernel : `numpy.ndarray`
    Simulated brighter-fatter kernel.

Definition at line 337 of file isrMock.py.

337 def makeBfKernel(self):
338 """Generate a simple Gaussian brighter-fatter kernel.
339
340 Returns
341 -------
342 kernel : `numpy.ndarray`
343 Simulated brighter-fatter kernel.
344 """
345 return self.bfKernel
346

◆ makeCrosstalkCoeff()

lsst.ip.isr.isrMock.IsrMock.makeCrosstalkCoeff ( self)
Generate the simulated crosstalk coefficients.

Returns
-------
coeffs : `numpy.ndarray`
    Simulated crosstalk coefficients.

Definition at line 358 of file isrMock.py.

358 def makeCrosstalkCoeff(self):
359 """Generate the simulated crosstalk coefficients.
360
361 Returns
362 -------
363 coeffs : `numpy.ndarray`
364 Simulated crosstalk coefficients.
365 """
366
367 return self.crosstalkCoeffs
368

◆ makeData()

lsst.ip.isr.isrMock.IsrMock.makeData ( self)
Generate simulated ISR data.

Currently, only the class defined crosstalk coefficient
matrix, brighter-fatter kernel, a constant unity transmission
curve, or a simple single-entry defect list can be generated.

Returns
-------
dataProduct :
    Simulated ISR data product.

Definition at line 309 of file isrMock.py.

309 def makeData(self):
310 """Generate simulated ISR data.
311
312 Currently, only the class defined crosstalk coefficient
313 matrix, brighter-fatter kernel, a constant unity transmission
314 curve, or a simple single-entry defect list can be generated.
315
316 Returns
317 -------
318 dataProduct :
319 Simulated ISR data product.
320 """
321 if sum(map(bool, [self.config.doBrighterFatter,
322 self.config.doDefects,
323 self.config.doTransmissionCurve,
324 self.config.doCrosstalkCoeffs])) != 1:
325 raise RuntimeError("Only one data product can be generated at a time.")
326 elif self.config.doBrighterFatter is True:
327 return self.makeBfKernel()
328 elif self.config.doDefects is True:
329 return self.makeDefectList()
330 elif self.config.doTransmissionCurve is True:
331 return self.makeTransmissionCurve()
332 elif self.config.doCrosstalkCoeffs is True:
333 return self.crosstalkCoeffs
334 else:
335 return None
336

◆ makeDefectList()

lsst.ip.isr.isrMock.IsrMock.makeDefectList ( self)
Generate a simple single-entry defect list.

Returns
-------
defectList : `lsst.meas.algorithms.Defects`
    Simulated defect list

Definition at line 347 of file isrMock.py.

347 def makeDefectList(self):
348 """Generate a simple single-entry defect list.
349
350 Returns
351 -------
352 defectList : `lsst.meas.algorithms.Defects`
353 Simulated defect list
354 """
355 return Defects([lsst.geom.Box2I(lsst.geom.Point2I(0, 0),
356 lsst.geom.Extent2I(40, 50))])
357
An integer coordinate rectangle.
Definition Box.h:55

◆ makeImage()

lsst.ip.isr.isrMock.IsrMock.makeImage ( self)
Generate a simulated ISR image.

Returns
-------
exposure : `lsst.afw.image.Exposure` or `dict`
    Simulated ISR image data.

Notes
-----
This method currently constructs a "raw" data image by:

* Generating a simulated sky with noise
* Adding a single Gaussian "star"
* Adding the fringe signal
* Multiplying the frame by the simulated flat
* Adding dark current (and noise)
* Adding a bias offset (and noise)
* Adding an overscan gradient parallel to the pixel y-axis
* Simulating crosstalk by adding a scaled version of each
  amplifier to each other amplifier.

The exposure with image data constructed this way is in one of
three formats.

* A single image, with overscan and prescan regions retained
* A single image, with overscan and prescan regions trimmed
* A `dict`, containing the amplifer data indexed by the
  amplifier name.

The nonlinearity, CTE, and brighter fatter are currently not
implemented.

Note that this method generates an image in the reverse
direction as the ISR processing, as the output image here has
had a series of instrument effects added to an idealized
exposure.

Definition at line 380 of file isrMock.py.

380 def makeImage(self):
381 """Generate a simulated ISR image.
382
383 Returns
384 -------
385 exposure : `lsst.afw.image.Exposure` or `dict`
386 Simulated ISR image data.
387
388 Notes
389 -----
390 This method currently constructs a "raw" data image by:
391
392 * Generating a simulated sky with noise
393 * Adding a single Gaussian "star"
394 * Adding the fringe signal
395 * Multiplying the frame by the simulated flat
396 * Adding dark current (and noise)
397 * Adding a bias offset (and noise)
398 * Adding an overscan gradient parallel to the pixel y-axis
399 * Simulating crosstalk by adding a scaled version of each
400 amplifier to each other amplifier.
401
402 The exposure with image data constructed this way is in one of
403 three formats.
404
405 * A single image, with overscan and prescan regions retained
406 * A single image, with overscan and prescan regions trimmed
407 * A `dict`, containing the amplifer data indexed by the
408 amplifier name.
409
410 The nonlinearity, CTE, and brighter fatter are currently not
411 implemented.
412
413 Note that this method generates an image in the reverse
414 direction as the ISR processing, as the output image here has
415 had a series of instrument effects added to an idealized
416 exposure.
417 """
418 exposure = self.getExposure()
419
420 for idx, amp in enumerate(exposure.getDetector()):
421 bbox = None
422 if self.config.isTrimmed is True:
423 bbox = amp.getBBox()
424 else:
425 bbox = amp.getRawDataBBox()
426
427 ampData = exposure.image[bbox]
428
429 if self.config.doAddSky is True:
430 self.amplifierAddNoise(ampData, self.config.skyLevel, np.sqrt(self.config.skyLevel))
431
432 if self.config.doAddSource is True:
433 for sourceAmp, sourceFlux, sourceX, sourceY in zip(self.config.sourceAmp,
434 self.config.sourceFlux,
435 self.config.sourceX,
436 self.config.sourceY):
437 if idx == sourceAmp:
438 self.amplifierAddSource(ampData, sourceFlux, sourceX, sourceY)
439
440 if self.config.doAddFringe is True:
441 self.amplifierAddFringe(amp, ampData, np.array(self.config.fringeScale),
442 x0=np.array(self.config.fringeX0),
443 y0=np.array(self.config.fringeY0))
444
445 if self.config.doAddFlat is True:
446 if ampData.getArray().sum() == 0.0:
447 self.amplifierAddNoise(ampData, 1.0, 0.0)
448 u0 = exposure.getDimensions().getX()
449 v0 = exposure.getDimensions().getY()
450 self.amplifierMultiplyFlat(amp, ampData, self.config.flatDrop, u0=u0, v0=v0)
451
452 if self.config.doAddDark is True:
453 self.amplifierAddNoise(ampData,
454 self.config.darkRate * self.config.darkTime / self.config.gain,
455 np.sqrt(self.config.darkRate
456 * self.config.darkTime / self.config.gain))
457
458 if self.config.doAddCrosstalk is True:
459 ctCalib = CrosstalkCalib()
460 for idxS, ampS in enumerate(exposure.getDetector()):
461 for idxT, ampT in enumerate(exposure.getDetector()):
462 ampDataT = exposure.image[ampT.getBBox()
463 if self.config.isTrimmed else ampT.getRawDataBBox()]
464 outAmp = ctCalib.extractAmp(exposure.getImage(), ampS, ampT,
465 isTrimmed=self.config.isTrimmed)
466 self.amplifierAddCT(outAmp, ampDataT, self.crosstalkCoeffs[idxS][idxT])
467
468 for amp in exposure.getDetector():
469 bbox = None
470 if self.config.isTrimmed is True:
471 bbox = amp.getBBox()
472 else:
473 bbox = amp.getRawDataBBox()
474
475 ampData = exposure.image[bbox]
476
477 if self.config.doAddBias is True:
478 self.amplifierAddNoise(ampData, self.config.biasLevel,
479 self.config.readNoise / self.config.gain)
480
481 if self.config.doAddOverscan is True:
482 oscanBBox = amp.getRawHorizontalOverscanBBox()
483 oscanData = exposure.image[oscanBBox]
484 self.amplifierAddNoise(oscanData, self.config.biasLevel,
485 self.config.readNoise / self.config.gain)
486
487 self.amplifierAddYGradient(ampData, -1.0 * self.config.overscanScale,
488 1.0 * self.config.overscanScale)
489 self.amplifierAddYGradient(oscanData, -1.0 * self.config.overscanScale,
490 1.0 * self.config.overscanScale)
491
492 if self.config.doGenerateAmpDict is True:
493 expDict = dict()
494 for amp in exposure.getDetector():
495 expDict[amp.getName()] = exposure
496 return expDict
497 else:
498 return exposure
499

◆ makeTransmissionCurve()

lsst.ip.isr.isrMock.IsrMock.makeTransmissionCurve ( self)
Generate a simulated flat transmission curve.

Returns
-------
transmission : `lsst.afw.image.TransmissionCurve`
    Simulated transmission curve.

Definition at line 369 of file isrMock.py.

369 def makeTransmissionCurve(self):
370 """Generate a simulated flat transmission curve.
371
372 Returns
373 -------
374 transmission : `lsst.afw.image.TransmissionCurve`
375 Simulated transmission curve.
376 """
377
378 return afwImage.TransmissionCurve.makeIdentity()
379

◆ run()

lsst.ip.isr.isrMock.IsrMock.run ( self)
Generate a mock ISR product, and return it.

Returns
-------
image : `lsst.afw.image.Exposure`
    Simulated ISR image with signals added.
dataProduct :
    Simulated ISR data products.
None :
    Returned if no valid configuration was found.

Raises
------
RuntimeError
    Raised if both doGenerateImage and doGenerateData are specified.

Definition at line 283 of file isrMock.py.

283 def run(self):
284 """Generate a mock ISR product, and return it.
285
286 Returns
287 -------
288 image : `lsst.afw.image.Exposure`
289 Simulated ISR image with signals added.
290 dataProduct :
291 Simulated ISR data products.
292 None :
293 Returned if no valid configuration was found.
294
295 Raises
296 ------
297 RuntimeError
298 Raised if both doGenerateImage and doGenerateData are specified.
299 """
300 if self.config.doGenerateImage and self.config.doGenerateData:
301 raise RuntimeError("Only one of doGenerateImage and doGenerateData may be specified.")
302 elif self.config.doGenerateImage:
303 return self.makeImage()
304 elif self.config.doGenerateData:
305 return self.makeData()
306 else:
307 return None
308

Member Data Documentation

◆ _DefaultName

str lsst.ip.isr.isrMock.IsrMock._DefaultName = "isrMock"
staticprotected

Definition at line 263 of file isrMock.py.

◆ bfKernel

lsst.ip.isr.isrMock.IsrMock.bfKernel

Definition at line 277 of file isrMock.py.

◆ ConfigClass

lsst.ip.isr.isrMock.IsrMock.ConfigClass = IsrMockConfig
static

Definition at line 262 of file isrMock.py.

◆ crosstalkCoeffs

lsst.ip.isr.isrMock.IsrMock.crosstalkCoeffs

Definition at line 268 of file isrMock.py.

◆ rng

lsst.ip.isr.isrMock.IsrMock.rng

Definition at line 267 of file isrMock.py.


The documentation for this class was generated from the following file: