28 from .
import SingleGaussianPsf
31 def plantSources(bbox, kwid, sky, coordList, addPoissonNoise=True):
32 """Make an exposure with stars (modelled as Gaussians)
34 @param bbox: parent bbox of exposure
35 @param kwid: kernel width (and height; kernel is square)
36 @param sky: amount of sky background (counts)
37 @param coordList: a list of [x, y, counts, sigma], where:
38 * x,y are relative to exposure origin
39 * counts is the integrated counts for the star
40 * sigma is the Gaussian sigma in pixels
41 @param addPoissonNoise: add Poisson noise to the exposure?
44 img = afwImage.ImageD(bbox)
46 for coord
in coordList:
47 x, y, counts, sigma = coord
51 psf = SingleGaussianPsf(kwid, kwid, sigma)
58 psfBox = thisPsfImg.getBBox()
60 if psfBox != thisPsfImg.getBBox():
61 thisPsfImg = thisPsfImg[psfBox, afwImage.PARENT]
62 imgSeg = img[psfBox, afwImage.PARENT]
64 meanSigma /= len(coordList)
70 np.random.seed(seed=1)
71 imgArr = img.getArray()
72 imgArr[:] = np.random.poisson(imgArr)
76 var = img.convertFloat()
78 mimg = afwImage.MaskedImageF(img.convertFloat(), mask, var)
82 psf = SingleGaussianPsf(kwid, kwid, meanSigma)
89 maxRadius=80.0, nRadii=30, perturb=0.05):
90 """Create a random TransmissionCurve with nontrivial spatial and
95 rng : numpy.random.RandomState
96 Random number generator.
98 Average minimum wavelength for generated TransmissionCurves (will be
100 maxWavelength : float
101 Average maximum wavelength for generated TransmissionCurves (will be
104 Number of samples in the wavelength dimension.
106 Average maximum radius for spatial variation (will be perturbed).
108 Number of samples in the radial dimension.
110 Fraction by which wavelength and radius bounds should be randomly
113 dWavelength = maxWavelength - minWavelength
115 def perturbed(x, s=perturb*dWavelength):
116 return x + 2.0*s*(rng.rand() - 0.5)
118 wavelengths = np.linspace(perturbed(minWavelength), perturbed(maxWavelength), nWavelengths)
119 radii = np.linspace(0.0, perturbed(maxRadius, perturb*maxRadius), nRadii)
120 throughput = np.zeros(wavelengths.shape + radii.shape, dtype=float)
123 peak0 = perturbed(0.9, 0.05)
124 start0 = perturbed(minWavelength + 0.25*dWavelength)
125 stop0 = perturbed(minWavelength + 0.75*dWavelength)
126 for i, r
in enumerate(radii):
127 mask = np.logical_and(wavelengths >= start0 + r, wavelengths <= stop0 + r)
128 throughput[mask, i] = peak0*(1.0 - r/1000.0)
129 return afwImage.TransmissionCurve.makeRadial(throughput, wavelengths, radii)