LSST Applications  21.0.0-172-gfb10e10a+18fedfabac,22.0.0+297cba6710,22.0.0+80564b0ff1,22.0.0+8d77f4f51a,22.0.0+a28f4c53b1,22.0.0+dcf3732eb2,22.0.1-1-g7d6de66+2a20fdde0d,22.0.1-1-g8e32f31+297cba6710,22.0.1-1-geca5380+7fa3b7d9b6,22.0.1-12-g44dc1dc+2a20fdde0d,22.0.1-15-g6a90155+515f58c32b,22.0.1-16-g9282f48+790f5f2caa,22.0.1-2-g92698f7+dcf3732eb2,22.0.1-2-ga9b0f51+7fa3b7d9b6,22.0.1-2-gd1925c9+bf4f0e694f,22.0.1-24-g1ad7a390+a9625a72a8,22.0.1-25-g5bf6245+3ad8ecd50b,22.0.1-25-gb120d7b+8b5510f75f,22.0.1-27-g97737f7+2a20fdde0d,22.0.1-32-gf62ce7b1+aa4237961e,22.0.1-4-g0b3f228+2a20fdde0d,22.0.1-4-g243d05b+871c1b8305,22.0.1-4-g3a563be+32dcf1063f,22.0.1-4-g44f2e3d+9e4ab0f4fa,22.0.1-42-gca6935d93+ba5e5ca3eb,22.0.1-5-g15c806e+85460ae5f3,22.0.1-5-g58711c4+611d128589,22.0.1-5-g75bb458+99c117b92f,22.0.1-6-g1c63a23+7fa3b7d9b6,22.0.1-6-g50866e6+84ff5a128b,22.0.1-6-g8d3140d+720564cf76,22.0.1-6-gd805d02+cc5644f571,22.0.1-8-ge5750ce+85460ae5f3,master-g6e05de7fdc+babf819c66,master-g99da0e417a+8d77f4f51a,w.2021.48
LSST Data Management Base Package
Public Member Functions | Public Attributes | Static Public Attributes | List of all members
lsst.ip.isr.isrMock.FringeMock Class Reference
Inheritance diagram for lsst.ip.isr.isrMock.FringeMock:
lsst.ip.isr.isrMock.MasterMock lsst.ip.isr.isrMock.IsrMock lsst.ip.isr.isrMock.UntrimmedFringeMock

Public Member Functions

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

Public Attributes

 rng
 
 crosstalkCoeffs
 
 bfKernel
 

Static Public Attributes

 ConfigClass = IsrMockConfig
 

Detailed Description

Simulated master fringe calibration.

Definition at line 865 of file isrMock.py.

Constructor & Destructor Documentation

◆ __init__()

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

Reimplemented from lsst.ip.isr.isrMock.MasterMock.

Reimplemented in lsst.ip.isr.isrMock.UntrimmedFringeMock.

Definition at line 868 of file isrMock.py.

868  def __init__(self, **kwargs):
869  super().__init__(**kwargs)
870  self.config.doAddFringe = True
871 
872 

Member Function Documentation

◆ amplifierAddCT()

def lsst.ip.isr.isrMock.IsrMock.amplifierAddCT (   self,
  ampDataSource,
  ampDataTarget,
  scale 
)
inherited
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 679 of file isrMock.py.

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

◆ amplifierAddFringe()

def lsst.ip.isr.isrMock.IsrMock.amplifierAddFringe (   self,
  amp,
  ampData,
  scale,
  x0 = 100,
  y0 = 0 
)
inherited
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 701 of file isrMock.py.

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

◆ amplifierAddNoise()

def lsst.ip.isr.isrMock.IsrMock.amplifierAddNoise (   self,
  ampData,
  mean,
  sigma 
)
inherited
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 621 of file isrMock.py.

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

◆ amplifierAddSource()

def lsst.ip.isr.isrMock.IsrMock.amplifierAddSource (   self,
  ampData,
  scale,
  x0,
  y0 
)
inherited
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 658 of file isrMock.py.

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

◆ amplifierAddYGradient()

def lsst.ip.isr.isrMock.IsrMock.amplifierAddYGradient (   self,
  ampData,
  start,
  end 
)
inherited
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 639 of file isrMock.py.

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

◆ amplifierMultiplyFlat()

def lsst.ip.isr.isrMock.IsrMock.amplifierMultiplyFlat (   self,
  amp,
  ampData,
  fracDrop,
  u0 = 100.0,
  v0 = 100.0 
)
inherited
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 732 of file isrMock.py.

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

◆ getCamera()

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

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

Definition at line 499 of file isrMock.py.

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

◆ getExposure()

def lsst.ip.isr.isrMock.IsrMock.getExposure (   self)
inherited
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 515 of file isrMock.py.

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

◆ getWcs()

def lsst.ip.isr.isrMock.IsrMock.getWcs (   self)
inherited
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 573 of file isrMock.py.

573  def getWcs(self):
574  """Construct a dummy WCS object.
575 
576  Taken from the deprecated ip_isr/examples/exampleUtils.py.
577 
578  This is not guaranteed, given the distortion and pixel scale
579  listed in the afwTestUtils camera definition.
580 
581  Returns
582  -------
583  wcs : `lsst.afw.geom.SkyWcs`
584  Test WCS transform.
585  """
586  return afwGeom.makeSkyWcs(crpix=lsst.geom.Point2D(0.0, 100.0),
587  crval=lsst.geom.SpherePoint(45.0, 25.0, lsst.geom.degrees),
588  cdMatrix=afwGeom.makeCdMatrix(scale=1.0*lsst.geom.degrees))
589 
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()

def lsst.ip.isr.isrMock.IsrMock.localCoordToExpCoord (   self,
  ampData,
  x,
  y 
)
inherited
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 590 of file isrMock.py.

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

◆ makeBfKernel()

def lsst.ip.isr.isrMock.IsrMock.makeBfKernel (   self)
inherited
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()

def lsst.ip.isr.isrMock.IsrMock.makeCrosstalkCoeff (   self)
inherited
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()

def lsst.ip.isr.isrMock.IsrMock.makeData (   self)
inherited
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()

def lsst.ip.isr.isrMock.IsrMock.makeDefectList (   self)
inherited
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()

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

◆ makeTransmissionCurve()

def lsst.ip.isr.isrMock.IsrMock.makeTransmissionCurve (   self)
inherited
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()

def lsst.ip.isr.isrMock.IsrMock.run (   self)
inherited
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 
def run(self, coaddExposures, bbox, wcs)
Definition: getTemplate.py:603

Member Data Documentation

◆ bfKernel

lsst.ip.isr.isrMock.IsrMock.bfKernel
inherited

Definition at line 277 of file isrMock.py.

◆ ConfigClass

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

Definition at line 262 of file isrMock.py.

◆ crosstalkCoeffs

lsst.ip.isr.isrMock.IsrMock.crosstalkCoeffs
inherited

Definition at line 268 of file isrMock.py.

◆ rng

lsst.ip.isr.isrMock.IsrMock.rng
inherited

Definition at line 267 of file isrMock.py.


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