

In [1]:

    
import galsim
import time
import numpy as np
import matplotlib.pyplot as plt

Test the accuracy of flux for truncated Sersic model in GalSim

Decide the parameters for Sersic component



In [37]:

    
flux = 1.00e6
nser = 4.0   
reff = 10.0   
ba   = 0.6 
pa   = 45.0

Decide the parameter for truncation



In [38]:

    
trunc = int(10.0 * reff)
method = 'auto'

Can also try to draw the object on an Image Object with the size we want



In [39]:

    
size = min([30.0 * reff, 600.0])
imgTemp = galsim.image.Image(size, size)

Make a Simple PSF for convolution test



In [40]:

    
psfObj = galsim.Moffat(2.0, fwhm=2.0)
plt.figure(1, figsize=(5,5))
plt.imshow(numpy.arcsinh(psfObj.drawImage().array))









    Out[40]:





<matplotlib.image.AxesImage at 0x2b60fba96410>

First, test the un-truncated Sersic model



In [21]:

    
start_time = time.time()

serUntr = galsim.Sersic(nser, half_light_radius=reff)
serUntr = serUntr.withFlux(flux)
print "\n## Normal round Sersic model with Sersic %5.1f, and Re %6.2f pixels" % (nser, reff)
print "The Flux is : ", serUntr.getFlux()
serUntrImg0 = serUntr.drawImage(method=method).array
print "Shape of the output array : ", serUntrImg0.shape
print "Total Flux in the array   : ", np.sum(serUntrImg0)

# Apply the axis ratio 
serUntr = serUntr.shear(q=ba, beta=0.0*galsim.degrees)
print "\n## With axis ratio = %5.2f" % ba
print "The Flux is : ", serUntr.getFlux()
serUntrImg1 = serUntr.drawImage().array
print "Shape of the output array : ", serUntrImg1.shape
print "Total Flux in the array   : ", np.sum(serUntrImg1)

# Apply the position angle
serUntr = serUntr.rotate(pa*galsim.degrees)
print "\n## With position angle = %6.2f" % pa
print "The Flux is : ", serUntr.getFlux()
serUntrImg2 = serUntr.drawImage(method=method).array
print "Shape of the output array : ", serUntrImg2.shape
print "Total Flux in the array   : ", np.sum(serUntrImg2)

# Apply the PSF convolution 
serConvUntr = galsim.Convolve([serUntr, psfObj])
print "\n## Convolve with a Moffat PSF model"
print "The Flux is : ", serConvUntr.getFlux()
serUntrImg3 = serConvUntr.drawImage(method=method).array
print "Shape of the output array : ", serUntrImg3.shape
print "Total Flux in the array   : ", np.sum(serUntrImg3)

print "\n------- Time Spent: %s seconds -------" % (time.time() - start_time)









    



## Truncated round Sersic model with Sersic   3.0, and Re   2.50 pixels
The Flux is :  1000000.0
Shape of the output array :  (1068, 1068)
Total Flux in the array   :  996619.0

## With axis ratio =  0.60
The Flux is :  1000000.0
Shape of the output array :  (1778, 1778)
Total Flux in the array   :  998602.0

## With position angle =  45.00
The Flux is :  1000000.0
Shape of the output array :  (1778, 1778)
Total Flux in the array   :  998772.0

## Convolve with a Moffat PSF model
The Flux is :  1000000.0
Shape of the output array :  (242, 242)
Total Flux in the array   :  999340.0

------- Time Spent: 0.611943006516 seconds -------

Show the results:



In [23]:

    
fig, axes = plt.subplots(1, 4, figsize=(16.0, 4.0))
imgs = serUntrImg0, serUntrImg1, serUntrImg2, serUntrImg3
titles = "SerOrigin", "SerAxisRatio", "SerPosAngle", "SerConv"
for i in range(4):
        axes[i].imshow(numpy.arcsinh(imgs[i]))
        axes[i].set_title(titles[i])

Second, test the un-truncated Sersic model



In [24]:

    
start_time = time.time()

serTrun = galsim.Sersic(nser, half_light_radius=reff, trunc=trunc)
serTrun = serTrun.withFlux(flux)
print "\n## Truncated round Sersic model with Sersic %5.1f, and Re %6.2f pixels" % (nser, reff)
print "The Flux is : ", serTrun.getFlux()
serTrunImg0 = serTrun.drawImage(method=method).array
print "Shape of the output array : ", serTrunImg0.shape
print "Total Flux in the array   : ", np.sum(serTrunImg0)

# Apply the axis ratio 
serTrun = serTrun.shear(q=ba, beta=0.0*galsim.degrees)
print "\n## With axis ratio = %5.2f" % ba
print "The Flux is : ", serTrun.getFlux()
serTrunImg1 = serTrun.drawImage().array
print "Shape of the output array : ", serTrunImg1.shape
print "Total Flux in the array   : ", np.sum(serTrunImg1)

# Apply the position angle
serTrun = serTrun.rotate(pa*galsim.degrees)
print "\n## With position angle = %6.2f" % pa
print "The Flux is : ", serTrun.getFlux()
serTrunImg2 = serTrun.drawImage(method=method).array
print "Shape of the output array : ", serTrunImg2.shape
print "Total Flux in the array   : ", np.sum(serTrunImg2)

# Apply the PSF convolution 
serConvTrun = galsim.Convolve([serTrun, psfObj])
print "\n## Convolve with a Moffat PSF model"
print "The Flux is : ", serConvTrun.getFlux()
serTrunImg3 = serConvTrun.drawImage(method=method).array
print "Shape of the output array : ", serTrunImg3.shape
print "Total Flux in the array   : ", np.sum(serTrunImg3)

print "\n------- Time Spent: %s seconds -------" % (time.time() - start_time)









    



## Normal round Sersic model with Sersic   3.0, and Re   2.50 pixels
The Flux is :  1000000.0
Shape of the output array :  (662, 662)
Total Flux in the array   :  1e+06

## With axis ratio =  0.60
The Flux is :  1000000.0
Shape of the output array :  (1102, 1102)
Total Flux in the array   :  1e+06

## With position angle =  45.00
The Flux is :  1000000.0
Shape of the output array :  (1102, 1102)
Total Flux in the array   :  1e+06

## Convolve with a Moffat PSF model
The Flux is :  1000000.0
Shape of the output array :  (174, 174)
Total Flux in the array   :  999469.0

------- Time Spent: 94.6557378769 seconds -------



In [25]:

    
fig, axes = plt.subplots(1, 4, figsize=(16.0, 4.0))
imgs = serTrunImg0, serTrunImg1, serTrunImg2, serTrunImg3
titles = "SerOrigin", "SerAxisRatio", "SerPosAngle", "SerConv"
for i in range(4):
        axes[i].imshow(numpy.arcsinh(imgs[i]))
        axes[i].set_title(titles[i])

Untruncated Case, Ignore the middle steps



In [31]:

    
start_time = time.time()

serUntr = galsim.Sersic(nser, half_light_radius=reff)
serUntr = serUntr.withFlux(flux)
print "\n## Normal round Sersic model with Sersic %5.1f, and Re %6.2f pixels" % (nser, reff)
print "The Flux is : ", serUntr.getFlux()

# Apply the axis ratio 
print "\n## With axis ratio = %5.2f" % ba
serUntr = serUntr.shear(q=ba, beta=0.0*galsim.degrees)

# Apply the position angle
print "\n## With position angle = %6.2f" % pa
serUntr = serUntr.rotate(pa*galsim.degrees)

# Apply the PSF convolution 
serConvUntr = galsim.Convolve([serUntr, psfObj])
print "\n## Convolve with a Moffat PSF model"
print "The Flux is : ", serConvUntr.getFlux()
serUntrImg = serConvUntr.drawImage(method=method).array
print "Shape of the output array : ", serUntrImg.shape
print "Total Flux in the array   : ", np.sum(serUntrImg)

print "\n------- Time Spent: %s seconds -------" % (time.time() - start_time)

plt.figure(1, figsize=(5,5))
plt.imshow(numpy.arcsinh(serUntrImg))









    



## Normal round Sersic model with Sersic   4.0, and Re  10.00 pixels
The Flux is :  1000000.0

## With axis ratio =  0.60

## With position angle =  45.00

## Convolve with a Moffat PSF model
The Flux is :  1000000.0
Shape of the output array :  (1414, 1414)
Total Flux in the array   :  999070.0

------- Time Spent: 1.13336491585 seconds -------






    Out[31]:





<matplotlib.image.AxesImage at 0x2b60fbedad90>

Truncated case, Ignore the middle steps



In [32]:

    
start_time = time.time()

serTrun = galsim.Sersic(nser, half_light_radius=reff, trunc=trunc)
serTrun = serTrun.withFlux(flux)
print "\n## Truncated round Sersic model with Sersic %5.1f, and Re %6.2f pixels" % (nser, reff)
print "The Flux is : ", serTrun.getFlux()

# Apply the axis ratio 
print "\n## With axis ratio = %5.2f" % ba
serTrun = serTrun.shear(q=ba, beta=0.0*galsim.degrees)

# Apply the position angle
print "\n## With position angle = %6.2f" % pa
serTrun = serTrun.rotate(pa*galsim.degrees)

# Apply the PSF convolution 
serConvTrun = galsim.Convolve([serTrun, psfObj])
print "\n## Convolve with a Moffat PSF model"
print "The Flux is : ", serConvTrun.getFlux()
serTrunImg = serConvTrun.drawImage(method=method).array
print "Shape of the output array : ", serTrunImg.shape
print "Total Flux in the array   : ", np.sum(serTrunImg)

print "\n------- Time Spent: %s seconds -------" % (time.time() - start_time)

plt.figure(1, figsize=(5,5))
plt.imshow(numpy.arcsinh(serTrunImg))









    



## Truncated round Sersic model with Sersic   4.0, and Re  10.00 pixels
The Flux is :  1000000.0

## With axis ratio =  0.60

## With position angle =  45.00

## Convolve with a Moffat PSF model
The Flux is :  1000000.0
Shape of the output array :  (174, 174)
Total Flux in the array   :  999275.0

------- Time Spent: 0.0334300994873 seconds -------






    Out[32]:





<matplotlib.image.AxesImage at 0x2b60fbe8d2d0>

Untruncated, but draw on an ImageObj with the size we want



In [42]:

    
start_time = time.time()

serUntr = galsim.Sersic(nser, half_light_radius=reff)
serUntr = serUntr.withFlux(flux)
print "\n## Normal round Sersic model with Sersic %5.1f, and Re %6.2f pixels" % (nser, reff)
print "The Flux is : ", serUntr.getFlux()

# Apply the axis ratio 
print "\n## With axis ratio = %5.2f" % ba
serUntr = serUntr.shear(q=ba, beta=0.0*galsim.degrees)

# Apply the position angle
print "\n## With position angle = %6.2f" % pa
serUntr = serUntr.rotate(pa*galsim.degrees)

# Apply the PSF convolution 
serConvUntr = galsim.Convolve([serUntr, psfObj])
print "\n## Convolve with a Moffat PSF model"
print "The Flux is : ", serConvUntr.getFlux()

# Draw on an ImgObj with size=size
print "\n## Draw on an ImageObj with size=%5d pixels" % size
serUntrImg = serConvUntr.drawImage(imgTemp, method=method).array
print "Shape of the output array : ", serUntrImg.shape
print "Total Flux in the array   : ", np.sum(serUntrImg)

print "\n------- Time Spent: %s seconds -------" % (time.time() - start_time)

plt.figure(1, figsize=(5,5))
plt.imshow(numpy.arcsinh(serUntrImg))









    



## Normal round Sersic model with Sersic   4.0, and Re  10.00 pixels
The Flux is :  1000000.0

## With axis ratio =  0.60

## With position angle =  45.00

## Convolve with a Moffat PSF model
The Flux is :  1000000.0

## Draw on an ImageObj with size=  300 pixels
Shape of the output array :  (300, 300)
Total Flux in the array   :  933883.0

------- Time Spent: 0.246865987778 seconds -------






    Out[42]:





<matplotlib.image.AxesImage at 0x2b60fb931f90>



In [ ]: