Commit 64c00169 authored by Fang Yuedong's avatar Fang Yuedong
Browse files

add unittest script for detect effects

parent e0fa93a6
Loading
Loading
Loading
Loading
+1 −1
Original line number Diff line number Diff line
@@ -41,7 +41,7 @@ class Catalog_example(CatalogBase):
        config : dict
            configuration dictionary which is parsed from the input YAML file
        chip: ObservationSim.Instrument.Chip
            a ObservationSim.Instrument.Chip object, can be used to identify the band etc.
            a ObservationSim.Instrument.Chip instance, can be used to identify the band etc.
        **kwargs : dict
            other needed input parameters (in key-value pairs), please modify corresponding
            initialization call in "ObservationSim.py" as you need.
+204 −0
Original line number Diff line number Diff line
import unittest
import numpy as np
from ObservationSim.Instrument.Chip import Effects
import galsim
import matplotlib.pyplot as plt
import os,sys,math,copy
from numpy.random import Generator, PCG64
import warnings

warnings.filterwarnings("ignore", '.*Numba.*',)

width = 9216
height = 9232

if os.path.isdir('./output/'):
    pass
else:
    os.mkdir('./output/')

class DetTest(unittest.TestCase):
    def test_prnu(self):
        '''
        Unit test for PRNU. Expected result: a randomized GS image contains PRNU with sigma=0.01, mean=1.
        '''
        print('PRNU Test:')
        sigma = 0.01
        seed = 20210911
        prnuimg = Effects.PRNU_Img(width, height, sigma=sigma, seed=seed)
        meanval, stdval = np.mean(prnuimg.array), np.std(prnuimg.array)
        print('  Mean & STDDEV of PRNU image are %6.4f & %6.4f.' % (meanval, stdval))
        print('  PRNU Image Array:')
        print('  ',prnuimg.array)
        self.assertTrue(np.abs(meanval-1)<1e-6)
        self.assertTrue(np.abs(stdval-sigma)<0.002)
        print('\nUnit test for PRNU has been passed.')
        del prnuimg

    def test_dark(self):
        '''
        Test add dark current to image. Expected result: an image with dark current 3.4 e- and noise=1.844 e-.
        '''
        rng_poisson = galsim.BaseDeviate(20210911)
        dark_noise = galsim.DeviateNoise(galsim.PoissonDeviate(rng_poisson, 0.02*(150+0.5*40)))
        img = galsim.Image(200,200,dtype=np.float32, init_value=0)
        print('Initial Mean & STD = %6.3f & %6.3f' % (np.mean(img.array), np.std(img.array)))
        img.addNoise(dark_noise)
        meanval = np.mean(img.array)
        stdval = np.std(img.array)
        print('Dark added Mean & STD = %6.3f & %6.3f' % (meanval, stdval))
        self.assertTrue(np.abs(meanval-3.4)<0.05)
        self.assertTrue(np.abs(stdval-1.844)<0.02)
        print('\nUnit test for dark current has been passed.')
        del img

    def test_satu(self):
        '''
        Test saturation and bleeding. Expected result: an image with bleeding effect.
        '''
        img = galsim.Image(500,500,dtype=np.float32)
        star = galsim.Gaussian(flux=60e5,fwhm=3)
        img = star.drawImage(image=img,center=(150,200))
        # gal = galsim.Sersic(n=1, half_light_radius=3,flux=50e5)
        # img = gal.drawImage(image=img,center=(350,300))
        img.addNoise(galsim.GaussianNoise(sigma=7))
        # plt.imshow(img.array)
        # plt.show()
        filename1 = os.path.join('output','test_satu_initimg.fits')
        img.write(filename1)
        newimg = Effects.SaturBloom(img, fullwell=9e4)
        # plt.imshow(newimg.array)
        # plt.show()
        filename2 = os.path.join('output','test_satu_bleedimg.fits')
        newimg.write(filename2)
        del img,newimg, star

    def test_nonlinear(self):
        '''
        Test non-linear effect. Expected result: an image with non-linearity effect.
        '''
        imgarr = np.arange(1,9e4,4).reshape((150,150))
        img = galsim.Image(copy.deepcopy(imgarr))
        filename1 = os.path.join('output','test_nonlinear_initimg.fits')
        img.write(filename1)
        newimg = Effects.NonLinearity(img, beta1=5E-7, beta2=0)
        filename2 = os.path.join('output','test_nonlinear_finalimg.fits')
        newimg.write(filename2)
        plt.scatter(imgarr.flatten(), newimg.array.flatten(), s=2, alpha=0.5)
        plt.plot([-1e3,9e4],[-1e3,9e4],color='black', lw=1, ls='--')
        plt.xlabel('input (e-)')
        plt.ylabel('output (e-)')
        plt.savefig('./output/test_nonlinearity.png', dpi=200)
        plt.show()
        del img,newimg,imgarr

    def test_badpixel_HtrDtr(self):
        img = galsim.Image(500,500,init_value=1000)
        rgbadpix = Generator(PCG64(20210911))
        badfraction = 5E-5*(rgbadpix.random()*0.5+0.7)
        img = Effects.DefectivePixels(img, IfHotPix=True, IfDeadPix=True, fraction=badfraction, seed=20210911, biaslevel=0)
        img.write('./output/test_badpixel_HtrDtr.fits')
        del img
    def test_badpixel_HfsDtr(self):
        img = galsim.Image(500,500,init_value=1000)
        rgbadpix = Generator(PCG64(20210911))
        badfraction = 5E-5*(rgbadpix.random()*0.5+0.7)
        img = Effects.DefectivePixels(img, IfHotPix=False, IfDeadPix=True, fraction=badfraction, seed=20210911, biaslevel=0)
        img.write('./output/test_badpixel_HfsDtr.fits')
        del img
    def test_badpixel_HtrDfs(self):
        img = galsim.Image(500,500,init_value=1000)
        rgbadpix = Generator(PCG64(20210911))
        badfraction = 5E-5*(rgbadpix.random()*0.5+0.7)
        img = Effects.DefectivePixels(img, IfHotPix=True, IfDeadPix=False, fraction=badfraction, seed=20210911, biaslevel=0)
        img.write('./output/test_badpixel_HtrDfs.fits')
        del img
    def test_badpixel_HfsDfs(self):
        img = galsim.Image(500,500,init_value=1000)
        rgbadpix = Generator(PCG64(20210911))
        badfraction = 5E-5*(rgbadpix.random()*0.5+0.7)
        img = Effects.DefectivePixels(img, IfHotPix=False, IfDeadPix=False, fraction=badfraction, seed=20210911, biaslevel=0)
        img.write('./output/test_badpixel_HfsDfs.fits')
        del img

    def test_badlines(self):
        img = galsim.Image(500,500,init_value=-1000)
        img.addNoise(galsim.GaussianNoise(sigma=7))
        newimg = Effects.BadColumns(copy.deepcopy(img), seed=20210911)
        newimg.write('./output/test_badlines.fits')
        del newimg,img

    def test_cte(self):
        img = galsim.Image(200,200,init_value=1000)
        img.array[50,80] = 1e4
        img.array[150,150] = 3e4
        newimgcol = Effects.CTE_Effect(copy.deepcopy(img),direction='column')
        newimgrow = Effects.CTE_Effect(copy.deepcopy(img),direction='row')
        newimgcol.write('./output/test_ctecol.fits')
        newimgrow.write('./output/test_cterow.fits')
        del img,newimgcol,newimgrow

    def test_readnoise(self):
        img = galsim.Image(200,200,init_value=1000)
        seed = 20210911
        rng_readout = galsim.BaseDeviate(seed)
        readout_noise = galsim.GaussianNoise(rng=rng_readout, sigma=5)
        img.addNoise(readout_noise)
        img.write('./output/test_readnoise.fits')
        stdval = np.std(img.array)
        self.assertTrue(np.abs(stdval-5)<0.01*5)
        print('\nUnit test for readout noise has been passed.')
        del img

    def test_addbias(self):
        img = galsim.Image(200,200,init_value=0)
        img = Effects.AddBiasNonUniform16(img,bias_level=500, nsecy = 2, nsecx=8,seed=20210911)
        img.write('./output/test_addbias.fits')
        del img

    def test_apply16gains(self):
        img = galsim.Image(500,500,init_value=100)
        img = Effects.ApplyGainNonUniform16(img, gain=1.5, nsecy=2, nsecx=8, seed=202102)
        img.write("./output/test_apply16gains.fits")
        rightedge = int(500/8)*8
        print('gain=%6.2f' % 1.5)
        meanimg = np.mean(img.array[:,:rightedge])
        sigmaimg = np.std(img.array[:,:rightedge])
        print('mean, sigma = %6.2f, %6.2f' % (meanimg,sigmaimg))
        self.assertTrue(np.abs(meanimg-100/1.5)<1)
        self.assertTrue(np.abs(sigmaimg/meanimg-0.01)<0.001)
        print('\nUnit test for applying 16 channel gains has been passed.')
        del img


    def test_cosmicray(self):
        attachedSizes = np.loadtxt('../ObservationSim/Instrument/Chip/wfc-cr-attachpixel.dat')
        cr_map = Effects.produceCR_Map(
                xLen=500, yLen=500, exTime=150+0.5*40, 
                cr_pixelRatio=0.003*(1+0.5*40/150), 
                gain=1, attachedSizes=attachedSizes, seed=20210911)
        crimg = galsim.Image(cr_map)
        crimg.write('./output/test_cosmicray.fits')
        del cr_map,crimg

    def test_shutter(self):
        img = galsim.Image(5000,5000,init_value=1000)
        shuttimg = Effects.ShutterEffectArr(img, t_exp=150, t_shutter=1.3, dist_bearing=735, dt=1E-3)    # shutter effect normalized image for this chip
        img *= shuttimg
        img.write('./output/test_shutter.fits')
        del img

    def test_vignette(self):
        img = galsim.Image(2000,2000,init_value=1000)
        print(img.bounds)
        # # img.bounds = galsim.BoundsI(1, width, 1, height)
        img.setOrigin(10000,10000)
        flat_img = Effects.MakeFlatSmooth(img.bounds,20210911)
        flat_normal = flat_img / np.mean(flat_img.array)
        flat_normal.write('./output/test_vignette.fits')
        del flat_img,img,flat_normal



if __name__ == '__main__':
    unittest.main()
 No newline at end of file