Skip to content
GitLab
Commit fe70f366 authored by Fang Yuedong's avatar Fang Yuedong
Browse files

add unittest for focal plane

parent fcfee677
Show whitespace changes
Inline Side-by-side
ObservationSim/Instrument/Chip/Chip.py
View file @ fe70f366
......@@ -280,774 +280,774 @@ class Chip(FocalPlane):
noise = self.dark_noise * exptime + self.read_noise**2
return noise
def addEffects(self, config, img, chip_output, filt, ra_cen, dec_cen, img_rot, exptime=150., pointing_ID=0, timestamp_obs=1621915200, pointing_type='SCI', sky_map=None, post_flash_map=None, tel=None, logger=None):
# Set random seeds
SeedGainNonuni = int(config["random_seeds"]["seed_gainNonUniform"])
SeedBiasNonuni = int(config["random_seeds"]["seed_biasNonUniform"])
SeedRnNonuni = int(config["random_seeds"]["seed_rnNonUniform"])
SeedBadColumns = int(config["random_seeds"]["seed_badcolumns"])
SeedDefective = int(config["random_seeds"]["seed_defective"])
SeedCosmicRay = int(config["random_seeds"]["seed_CR"])
fullwell = int(self.full_well)
if config["ins_effects"]["add_hotpixels"] == True:
BoolHotPix = True
else:
BoolHotPix = False
if config["ins_effects"]["add_deadpixels"] == True:
BoolDeadPix = True
else:
BoolDeadPix = False
self.logger = logger
# Get Poisson noise generator
rng_poisson, poisson_noise = chip_utils.get_poisson(
seed=int(config["random_seeds"]["seed_poisson"]) + pointing_ID*30 + self.chipID, sky_level=0.)
# Add sky background
if config["ins_effects"]["add_back"] == True:
img, sky_map = chip_utils.add_sky_background(
img=img, filt=filt, exptime=exptime, sky_map=sky_map, tel=tel)
del sky_map
# Apply flat-field large scale structure for one chip
if config["ins_effects"]["flat_fielding"] == True:
chip_utils.log_info(
msg=" Creating and applying Flat-Fielding", logger=self.logger)
chip_utils.log_info(msg=str(img.bounds), logger=self.logger)
flat_img, flat_normal = chip_utils.get_flat(
img=img, seed=int(config["random_seeds"]["seed_flat"]))
if self.survey_type == "photometric":
img *= flat_normal
del flat_normal
if config["output_setting"]["flat_output"] == False:
del flat_img
if post_flash_map is not None:
img = img + post_flash_map
# Apply Shutter-effect for one chip
if config["ins_effects"]["shutter_effect"] == True:
chip_utils.log_info(
msg=" Apply shutter effect", logger=self.logger)
# shutter effect normalized image for this chip
shuttimg = effects.ShutterEffectArr(
img, t_shutter=1.3, dist_bearing=735, dt=1E-3)
if self.survey_type == "photometric":
img *= shuttimg
# output 16-bit shutter effect image with pixel value <=65535
if config["output_setting"]["shutter_output"] == True:
shutt_gsimg = galsim.ImageUS(shuttimg*6E4)
shutt_gsimg.write("%s/ShutterEffect_%s_1.fits" %
(chip_output.subdir, self.chipID))
del shutt_gsimg
del shuttimg
# # Add Poisson noise to the resulting images
# # (NOTE): this can only applied to the slitless image
# # since it dose not use photon shooting to draw stamps
# if self.survey_type == "spectroscopic":
# img.addNoise(poisson_noise)
# Add cosmic-rays
if config["ins_effects"]["cosmic_ray"] == True and pointing_type == 'SCI':
chip_utils.log_info(msg=" Adding Cosmic-Ray", logger=self.logger)
img, crmap_gsimg, cr_event_num = chip_utils.add_cosmic_rays(img=img, chip=self, exptime=exptime,
seed=SeedCosmicRay+pointing_ID*30+self.chipID)
chip_utils.outputCal(
chip=self,
img=crmap_gsimg,
ra_cen=ra_cen,
dec_cen=dec_cen,
img_rot=img_rot,
im_type='CRS',
pointing_ID=pointing_ID,
output_dir=chip_output.subdir,
exptime=exptime,
project_cycle=config["project_cycle"],
run_counter=config["run_counter"],
timestamp=timestamp_obs)
del crmap_gsimg
# Apply PRNU effect and output PRNU flat file:
if config["ins_effects"]["prnu_effect"] == True:
chip_utils.log_info(
msg=" Applying PRNU effect", logger=self.logger)
img, prnu_img = chip_utils.add_PRNU(img=img, chip=self,
seed=int(config["random_seeds"]["seed_prnu"]+self.chipID))
if config["output_setting"]["prnu_output"] == True:
prnu_img.write("%s/FlatImg_PRNU_%s.fits" %
(chip_output.subdir, self.chipID))
if config["output_setting"]["flat_output"] == False:
del prnu_img
# # Add dark current
# def addEffects(self, config, img, chip_output, filt, ra_cen, dec_cen, img_rot, exptime=150., pointing_ID=0, timestamp_obs=1621915200, pointing_type='SCI', sky_map=None, post_flash_map=None, tel=None, logger=None):
# # Set random seeds
# SeedGainNonuni = int(config["random_seeds"]["seed_gainNonUniform"])
# SeedBiasNonuni = int(config["random_seeds"]["seed_biasNonUniform"])
# SeedRnNonuni = int(config["random_seeds"]["seed_rnNonUniform"])
# SeedBadColumns = int(config["random_seeds"]["seed_badcolumns"])
# SeedDefective = int(config["random_seeds"]["seed_defective"])
# SeedCosmicRay = int(config["random_seeds"]["seed_CR"])
# fullwell = int(self.full_well)
# if config["ins_effects"]["add_hotpixels"] == True:
# BoolHotPix = True
# else:
# BoolHotPix = False
# if config["ins_effects"]["add_deadpixels"] == True:
# BoolDeadPix = True
# else:
# BoolDeadPix = False
# self.logger = logger
# # Get Poisson noise generator
# rng_poisson, poisson_noise = chip_utils.get_poisson(
# seed=int(config["random_seeds"]["seed_poisson"]) + pointing_ID*30 + self.chipID, sky_level=0.)
# # Add sky background
# if config["ins_effects"]["add_back"] == True:
# img, sky_map = chip_utils.add_sky_background(
# img=img, filt=filt, exptime=exptime, sky_map=sky_map, tel=tel)
# del sky_map
# # Apply flat-field large scale structure for one chip
# if config["ins_effects"]["flat_fielding"] == True:
# chip_utils.log_info(
# msg=" Creating and applying Flat-Fielding", logger=self.logger)
# chip_utils.log_info(msg=str(img.bounds), logger=self.logger)
# flat_img, flat_normal = chip_utils.get_flat(
# img=img, seed=int(config["random_seeds"]["seed_flat"]))
# if self.survey_type == "photometric":
# img *= flat_normal
# del flat_normal
# if config["output_setting"]["flat_output"] == False:
# del flat_img
# if post_flash_map is not None:
# img = img + post_flash_map
# # Apply Shutter-effect for one chip
# if config["ins_effects"]["shutter_effect"] == True:
# chip_utils.log_info(
# msg=" Apply shutter effect", logger=self.logger)
# # shutter effect normalized image for this chip
# shuttimg = effects.ShutterEffectArr(
# img, t_shutter=1.3, dist_bearing=735, dt=1E-3)
# if self.survey_type == "photometric":
# img *= shuttimg
# # output 16-bit shutter effect image with pixel value <=65535
# if config["output_setting"]["shutter_output"] == True:
# shutt_gsimg = galsim.ImageUS(shuttimg*6E4)
# shutt_gsimg.write("%s/ShutterEffect_%s_1.fits" %
# (chip_output.subdir, self.chipID))
# del shutt_gsimg
# del shuttimg
# # # Add Poisson noise to the resulting images
# # # (NOTE): this can only applied to the slitless image
# # # since it dose not use photon shooting to draw stamps
# # if self.survey_type == "spectroscopic":
# # img.addNoise(poisson_noise)
# # Add cosmic-rays
# if config["ins_effects"]["cosmic_ray"] == True and pointing_type == 'SCI':
# chip_utils.log_info(msg=" Adding Cosmic-Ray", logger=self.logger)
# img, crmap_gsimg, cr_event_num = chip_utils.add_cosmic_rays(img=img, chip=self, exptime=exptime,
# seed=SeedCosmicRay+pointing_ID*30+self.chipID)
# chip_utils.outputCal(
# chip=self,
# img=crmap_gsimg,
# ra_cen=ra_cen,
# dec_cen=dec_cen,
# img_rot=img_rot,
# im_type='CRS',
# pointing_ID=pointing_ID,
# output_dir=chip_output.subdir,
# exptime=exptime,
# project_cycle=config["project_cycle"],
# run_counter=config["run_counter"],
# timestamp=timestamp_obs)
# del crmap_gsimg
# # Apply PRNU effect and output PRNU flat file:
# if config["ins_effects"]["prnu_effect"] == True:
# chip_utils.log_info(
# msg=" Applying PRNU effect", logger=self.logger)
# img, prnu_img = chip_utils.add_PRNU(img=img, chip=self,
# seed=int(config["random_seeds"]["seed_prnu"]+self.chipID))
# if config["output_setting"]["prnu_output"] == True:
# prnu_img.write("%s/FlatImg_PRNU_%s.fits" %
# (chip_output.subdir, self.chipID))
# if config["output_setting"]["flat_output"] == False:
# del prnu_img
# # # Add dark current
# # if config["ins_effects"]["add_dark"] == True:
# # dark_noise = galsim.DeviateNoise(galsim.PoissonDeviate(rng_poisson, self.dark_noise*(exptime+0.5*self.readout_time)))
# # img.addNoise(dark_noise)
# # Add dark current & Poisson noise
# InputDark = False
# if config["ins_effects"]["add_dark"] == True:
# dark_noise = galsim.DeviateNoise(galsim.PoissonDeviate(rng_poisson, self.dark_noise*(exptime+0.5*self.readout_time)))
# img.addNoise(dark_noise)
# Add dark current & Poisson noise
InputDark = False
if config["ins_effects"]["add_dark"] == True:
if InputDark:
img = chip_utils.add_inputdark(
img=img, chip=self, exptime=exptime)
else:
img, _ = chip_utils.add_poisson(
img=img, chip=self, exptime=exptime, poisson_noise=poisson_noise)
else:
img, _ = chip_utils.add_poisson(
img=img, chip=self, exptime=exptime, poisson_noise=poisson_noise, dark_noise=0.)
# Add diffusion & brighter-fatter effects
if config["ins_effects"]["bright_fatter"] == True:
img = chip_utils.add_brighter_fatter(img=img)
# Add Hot Pixels or/and Dead Pixels
rgbadpix = Generator(PCG64(int(SeedDefective+self.chipID)))
badfraction = 5E-5*(rgbadpix.random()*0.5+0.7)
img = effects.DefectivePixels(img, IfHotPix=BoolHotPix, IfDeadPix=BoolDeadPix,
fraction=badfraction, seed=SeedDefective+self.chipID, biaslevel=0)
# Apply Bad lines
if config["ins_effects"]["add_badcolumns"] == True:
img = effects.BadColumns(
img, seed=SeedBadColumns, chipid=self.chipID, logger=self.logger)
# Apply Nonlinearity on the chip image
if config["ins_effects"]["non_linear"] == True:
chip_utils.log_info(
msg=" Applying Non-Linearity on the chip image", logger=self.logger)
img = effects.NonLinearity(GSImage=img, beta1=5.e-7, beta2=0)
# Apply CCD Saturation & Blooming
if config["ins_effects"]["saturbloom"] == True:
chip_utils.log_info(
msg=" Applying CCD Saturation & Blooming", logger=self.logger)
img = effects.SaturBloom(
GSImage=img, nsect_x=1, nsect_y=1, fullwell=fullwell)
# Apply CTE Effect
# if InputDark:
# img = chip_utils.add_inputdark(
# img=img, chip=self, exptime=exptime)
# else:
# img, _ = chip_utils.add_poisson(
# img=img, chip=self, exptime=exptime, poisson_noise=poisson_noise)
# else:
# img, _ = chip_utils.add_poisson(
# img=img, chip=self, exptime=exptime, poisson_noise=poisson_noise, dark_noise=0.)
# # Add diffusion & brighter-fatter effects
# if config["ins_effects"]["bright_fatter"] == True:
# img = chip_utils.add_brighter_fatter(img=img)
# # Add Hot Pixels or/and Dead Pixels
# rgbadpix = Generator(PCG64(int(SeedDefective+self.chipID)))
# badfraction = 5E-5*(rgbadpix.random()*0.5+0.7)
# img = effects.DefectivePixels(img, IfHotPix=BoolHotPix, IfDeadPix=BoolDeadPix,
# fraction=badfraction, seed=SeedDefective+self.chipID, biaslevel=0)
# # Apply Bad lines
# if config["ins_effects"]["add_badcolumns"] == True:
# img = effects.BadColumns(
# img, seed=SeedBadColumns, chipid=self.chipID, logger=self.logger)
# # Apply Nonlinearity on the chip image
# if config["ins_effects"]["non_linear"] == True:
# chip_utils.log_info(
# msg=" Applying Non-Linearity on the chip image", logger=self.logger)
# img = effects.NonLinearity(GSImage=img, beta1=5.e-7, beta2=0)
# # Apply CCD Saturation & Blooming
# if config["ins_effects"]["saturbloom"] == True:
# chip_utils.log_info(
# msg=" Applying CCD Saturation & Blooming", logger=self.logger)
# img = effects.SaturBloom(
# GSImage=img, nsect_x=1, nsect_y=1, fullwell=fullwell)
# # Apply CTE Effect
# # if config["ins_effects"]["cte_trail"] == True:
# # chip_utils.log_info(msg=" Apply CTE Effect", logger=self.logger)
# # img = effects.CTE_Effect(GSImage=img, threshold=27)
# pre1 = self.prescan_x # 27
# over1 = self.overscan_x # 71
# pre2 = self.prescan_y # 0 #4
# over2 = self.overscan_y # 84 #80
# if config["ins_effects"]["cte_trail"] == True:
# chip_utils.log_info(msg=" Apply CTE Effect", logger=self.logger)
# img = effects.CTE_Effect(GSImage=img, threshold=27)
pre1 = self.prescan_x # 27
over1 = self.overscan_x # 71
pre2 = self.prescan_y # 0 #4
over2 = self.overscan_y # 84 #80
if config["ins_effects"]["cte_trail"] == True:
chip_utils.log_info(msg=" Apply CTE Effect", logger=self.logger)
# img = effects.CTE_Effect(GSImage=img, threshold=27)
# CTI_modeling
# 2*8 -> 1*16 img-layout
img = chip_utils.formatOutput(GSImage=img)
self.nsecy = 1
self.nsecx = 16
img_arr = img.array
ny, nx = img_arr.shape
dx = int(nx/self.nsecx)
dy = int(ny/self.nsecy)
newimg = galsim.Image(nx, int(ny+over2), init_value=0)
for ichannel in range(16):
print('\n***add CTI effects: pointing-{:} chip-{:} channel-{:}***'.format(
pointing_ID, self.chipID, ichannel+1))
# nx,ny,noverscan,nsp,nmax = 4608,4616,84,3,10
noverscan, nsp, nmax = over2, 3, 10
beta, w, c = 0.478, 84700, 0
t = np.array([0.74, 7.7, 37], dtype=np.float32)
rho_trap = np.array([0.6, 1.6, 1.4], dtype=np.float32)
trap_seeds = np.array(
[0, 1000, 10000], dtype=np.int32) + ichannel + self.chipID*16
release_seed = 50 + ichannel + pointing_ID*30 + self.chipID*16
newimg.array[:, 0+ichannel*dx:dx+ichannel*dx] = CTI_sim(
img_arr[:, 0+ichannel*dx:dx+ichannel*dx], dx, dy, noverscan, nsp, nmax, beta, w, c, t, rho_trap, trap_seeds, release_seed)
newimg.wcs = img.wcs
del img
img = newimg
# 1*16 -> 2*8 img-layout
img = chip_utils.formatRevert(GSImage=img)
self.nsecy = 2
self.nsecx = 8
# prescan & overscan
if config["ins_effects"]["add_prescan"] == True:
chip_utils.log_info(
msg=" Apply pre/over-scan", logger=self.logger)
if config["ins_effects"]["cte_trail"] == False:
img = chip_utils.AddPreScan(
GSImage=img, pre1=pre1, pre2=pre2, over1=over1, over2=over2)
if config["ins_effects"]["cte_trail"] == True:
img = chip_utils.AddPreScan(
GSImage=img, pre1=pre1, pre2=pre2, over1=over1, over2=0)
# 1*16 output
if config["ins_effects"]["format_output"] == True:
chip_utils.log_info(msg=" Apply 1*16 format", logger=self.logger)
img = chip_utils.formatOutput(GSImage=img)
self.nsecy = 1
self.nsecx = 16
# Add Bias level
if config["ins_effects"]["add_bias"] == True:
chip_utils.log_info(
msg=" Adding Bias level and 16-channel non-uniformity", logger=self.logger)
if config["ins_effects"]["bias_16channel"] == True:
img = effects.AddBiasNonUniform16(img,
bias_level=float(
self.bias_level),
nsecy=self.nsecy, nsecx=self.nsecx,
seed=SeedBiasNonuni+self.chipID,
logger=self.logger)
elif config["ins_effects"]["bias_16channel"] == False:
img += self.bias_level
# Add Read-out Noise
if config["ins_effects"]["add_readout"] == True:
seed = int(config["random_seeds"]["seed_readout"]
) + pointing_ID*30 + self.chipID
rng_readout = galsim.BaseDeviate(seed)
readout_noise = galsim.GaussianNoise(
rng=rng_readout, sigma=self.read_noise)
img.addNoise(readout_noise)
# Apply Gain & Quantization
chip_utils.log_info(
msg=" Applying Gain (and 16 channel non-uniformity) & Quantization", logger=self.logger)
if config["ins_effects"]["gain_16channel"] == True:
img, self.gain_channel = effects.ApplyGainNonUniform16(
img, gain=self.gain,
nsecy=self.nsecy, nsecx=self.nsecx,
seed=SeedGainNonuni+self.chipID,
logger=self.logger)
elif config["ins_effects"]["gain_16channel"] == False:
img /= self.gain
img.array[img.array > 65535] = 65535
img.replaceNegative(replace_value=0)
img.quantize()
######################################################################################
# Output images for calibration pointing
######################################################################################
# Bias output
if config["ins_effects"]["add_bias"] == True and config["output_setting"]["bias_output"] == True and pointing_type == 'CAL':
if self.logger is not None:
self.logger.info(" Output N frame Bias files")
else:
print(" Output N frame Bias files", flush=True)
NBias = int(config["output_setting"]["NBias"])
for i in range(NBias):
# BiasCombImg, BiasTag = effects.MakeBiasNcomb(
# self.npix_x, self.npix_y,
# bias_level=float(self.bias_level),
# ncombine=1, read_noise=self.read_noise, gain=1,
# # img = effects.CTE_Effect(GSImage=img, threshold=27)
# # CTI_modeling
# # 2*8 -> 1*16 img-layout
# img = chip_utils.formatOutput(GSImage=img)
# self.nsecy = 1
# self.nsecx = 16
# img_arr = img.array
# ny, nx = img_arr.shape
# dx = int(nx/self.nsecx)
# dy = int(ny/self.nsecy)
# newimg = galsim.Image(nx, int(ny+over2), init_value=0)
# for ichannel in range(16):
# print('\n***add CTI effects: pointing-{:} chip-{:} channel-{:}***'.format(
# pointing_ID, self.chipID, ichannel+1))
# # nx,ny,noverscan,nsp,nmax = 4608,4616,84,3,10
# noverscan, nsp, nmax = over2, 3, 10
# beta, w, c = 0.478, 84700, 0
# t = np.array([0.74, 7.7, 37], dtype=np.float32)
# rho_trap = np.array([0.6, 1.6, 1.4], dtype=np.float32)
# trap_seeds = np.array(
# [0, 1000, 10000], dtype=np.int32) + ichannel + self.chipID*16
# release_seed = 50 + ichannel + pointing_ID*30 + self.chipID*16
# newimg.array[:, 0+ichannel*dx:dx+ichannel*dx] = CTI_sim(
# img_arr[:, 0+ichannel*dx:dx+ichannel*dx], dx, dy, noverscan, nsp, nmax, beta, w, c, t, rho_trap, trap_seeds, release_seed)
# newimg.wcs = img.wcs
# del img
# img = newimg
# # 1*16 -> 2*8 img-layout
# img = chip_utils.formatRevert(GSImage=img)
# self.nsecy = 2
# self.nsecx = 8
# # prescan & overscan
# if config["ins_effects"]["add_prescan"] == True:
# chip_utils.log_info(
# msg=" Apply pre/over-scan", logger=self.logger)
# if config["ins_effects"]["cte_trail"] == False:
# img = chip_utils.AddPreScan(
# GSImage=img, pre1=pre1, pre2=pre2, over1=over1, over2=over2)
# if config["ins_effects"]["cte_trail"] == True:
# img = chip_utils.AddPreScan(
# GSImage=img, pre1=pre1, pre2=pre2, over1=over1, over2=0)
# # 1*16 output
# if config["ins_effects"]["format_output"] == True:
# chip_utils.log_info(msg=" Apply 1*16 format", logger=self.logger)
# img = chip_utils.formatOutput(GSImage=img)
# self.nsecy = 1
# self.nsecx = 16
# # Add Bias level
# if config["ins_effects"]["add_bias"] == True:
# chip_utils.log_info(
# msg=" Adding Bias level and 16-channel non-uniformity", logger=self.logger)
# if config["ins_effects"]["bias_16channel"] == True:
# img = effects.AddBiasNonUniform16(img,
# bias_level=float(
# self.bias_level),
# nsecy=self.nsecy, nsecx=self.nsecx,
# seed=SeedBiasNonuni+self.chipID,
# logger=self.logger)
BiasCombImg = galsim.Image(
self.npix_x, self.npix_y, init_value=0)
if config["ins_effects"]["add_bias"] == True:
biaslevel = self.bias_level
overscan = biaslevel-2
elif config["ins_effects"]["add_bias"] == False:
biaslevel = 0
overscan = 0
# Readout noise for Biases is not generated with random seeds. So readout noise for bias images can't be reproduced.
if config["ins_effects"]["cosmic_ray"] == True:
if config["ins_effects"]["cray_differ"] == True:
cr_map, cr_event_num = effects.produceCR_Map(
xLen=self.npix_x, yLen=self.npix_y,
exTime=0.01,
cr_pixelRatio=0.003 *
(0.01+0.5*self.readout_time)/150.,
gain=self.gain,
attachedSizes=self.attachedSizes,
seed=SeedCosmicRay+pointing_ID*30+self.chipID+1)
# seed: obj-imaging:+0; bias:+1; dark:+2; flat:+3;
BiasCombImg += cr_map
del cr_map
# Apply Bad lines
if config["ins_effects"]["add_badcolumns"] == True:
BiasCombImg = effects.BadColumns(
BiasCombImg-float(self.bias_level)+5, seed=SeedBadColumns, chipid=self.chipID, logger=self.logger) + float(self.bias_level)-5
# Non-Linearity for Bias
if config["ins_effects"]["non_linear"] == True:
if self.logger is not None:
self.logger.info(
" Applying Non-Linearity on the Bias image")
else:
print(
" Applying Non-Linearity on the Bias image", flush=True)
BiasCombImg = effects.NonLinearity(
GSImage=BiasCombImg, beta1=5.e-7, beta2=0)
# START
pre1 = self.prescan_x # 27
over1 = self.overscan_x # 71
pre2 = self.prescan_y # 0 #4
over2 = self.overscan_y # 84 #80
# prescan & overscan
if config["ins_effects"]["add_prescan"] == True:
chip_utils.log_info(
msg=" Apply pre/over-scan", logger=self.logger)
BiasCombImg = chip_utils.AddPreScan(
GSImage=BiasCombImg, pre1=pre1, pre2=pre2, over1=over1, over2=over2)
# 1*16 output
if config["ins_effects"]["format_output"] == True:
chip_utils.log_info(
msg=" Apply 1*16 format", logger=self.logger)
BiasCombImg = chip_utils.formatOutput(GSImage=BiasCombImg)
self.nsecy = 1
self.nsecx = 16
# END
# Add Bias level
if config["ins_effects"]["add_bias"] == True:
if self.logger is not None:
self.logger.info(
" Adding Bias level and 16-channel non-uniformity")
else:
print(" Adding Bias level and 16-channel non-uniformity")
BiasCombImg = effects.AddBiasNonUniform16(BiasCombImg,
bias_level=biaslevel,
nsecy=self.nsecy, nsecx=self.nsecx,
seed=SeedBiasNonuni+self.chipID,
logger=self.logger)
rng = galsim.UniformDeviate()
ncombine = 1
NoiseBias = galsim.GaussianNoise(
rng=rng, sigma=self.read_noise*ncombine**0.5)
BiasCombImg.addNoise(NoiseBias)
BiasCombImg, self.gain_channel = effects.ApplyGainNonUniform16(BiasCombImg, gain=self.gain,
nsecy=self.nsecy, nsecx=self.nsecx,
seed=SeedGainNonuni+self.chipID,
logger=self.logger)
# BiasCombImg = effects.AddOverscan(
# BiasCombImg,
# overscan=float(config["ins_effects"]["bias_level"])-2, gain=self.gain,
# widthl=27, widthr=27, widtht=8, widthb=8)
BiasCombImg.replaceNegative(replace_value=0)
BiasCombImg.quantize()
BiasCombImg = galsim.ImageUS(BiasCombImg)
timestamp_obs += 10 * 60
chip_utils.outputCal(
chip=self,
img=BiasCombImg,
ra_cen=ra_cen,
dec_cen=dec_cen,
img_rot=img_rot,
im_type='BIAS',
pointing_ID=pointing_ID,
output_dir=chip_output.subdir,
exptime=0.0,
project_cycle=config["project_cycle"],
run_counter=config["run_counter"],
timestamp=timestamp_obs)
del BiasCombImg
# Export combined (ncombine, Vignetting + PRNU) & single vignetting flat-field file
if config["ins_effects"]["flat_fielding"] == True and config["output_setting"]["flat_output"] == True and pointing_type == 'CAL':
if self.logger is not None:
self.logger.info(" Output N frame Flat-Field files")
else:
print(" Output N frame Flat-Field files", flush=True)
NFlat = int(config["output_setting"]["NFlat"])
if config["ins_effects"]["add_bias"] == True:
biaslevel = self.bias_level
overscan = biaslevel-2
elif config["ins_effects"]["add_bias"] == False:
biaslevel = 0
overscan = 0
darklevel = self.dark_noise * \
(self.flat_exptime+0.5*self.readout_time)
for i in range(NFlat):
FlatSingle = flat_img * prnu_img + darklevel
FlatCombImg, FlatTag = effects.MakeFlatNcomb(
flat_single_image=FlatSingle,
ncombine=1,
read_noise=self.read_noise,
gain=1,
overscan=overscan,
biaslevel=0,
seed_bias=SeedDefective+self.chipID,
logger=self.logger
)
if config["ins_effects"]["cosmic_ray"] == True:
if config["ins_effects"]["cray_differ"] == True:
cr_map, cr_event_num = effects.produceCR_Map(
xLen=self.npix_x, yLen=self.npix_y,
exTime=self.flat_exptime+0.5*self.readout_time,
cr_pixelRatio=0.003 *
(self.flat_exptime+0.5*self.readout_time)/150.,
gain=self.gain,
attachedSizes=self.attachedSizes,
seed=SeedCosmicRay+pointing_ID*30+self.chipID+3)
# seed: obj-imaging:+0; bias:+1; dark:+2; flat:+3;
FlatCombImg += cr_map
del cr_map
# Add Hot Pixels or/and Dead Pixels
rgbadpix = Generator(PCG64(int(SeedDefective+self.chipID)))
badfraction = 5E-5*(rgbadpix.random()*0.5+0.7)
FlatCombImg = effects.DefectivePixels(
FlatCombImg, IfHotPix=BoolHotPix, IfDeadPix=BoolDeadPix, fraction=badfraction, seed=SeedDefective+self.chipID, biaslevel=0)
# Apply Bad lines
if config["ins_effects"]["add_badcolumns"] == True:
FlatCombImg = effects.BadColumns(
FlatCombImg, seed=SeedBadColumns, chipid=self.chipID, logger=self.logger)
if config["ins_effects"]["non_linear"] == True:
if self.logger is not None:
self.logger.info(
" Applying Non-Linearity on the Flat image")
else:
print(
" Applying Non-Linearity on the Flat image", flush=True)
FlatCombImg = effects.NonLinearity(
GSImage=FlatCombImg, beta1=5.e-7, beta2=0)
# elif config["ins_effects"]["bias_16channel"] == False:
# img += self.bias_level
# # Add Read-out Noise
# if config["ins_effects"]["add_readout"] == True:
# seed = int(config["random_seeds"]["seed_readout"]
# ) + pointing_ID*30 + self.chipID
# rng_readout = galsim.BaseDeviate(seed)
# readout_noise = galsim.GaussianNoise(
# rng=rng_readout, sigma=self.read_noise)
# img.addNoise(readout_noise)
# # Apply Gain & Quantization
# chip_utils.log_info(
# msg=" Applying Gain (and 16 channel non-uniformity) & Quantization", logger=self.logger)
# if config["ins_effects"]["gain_16channel"] == True:
# img, self.gain_channel = effects.ApplyGainNonUniform16(
# img, gain=self.gain,
# nsecy=self.nsecy, nsecx=self.nsecx,
# seed=SeedGainNonuni+self.chipID,
# logger=self.logger)
# elif config["ins_effects"]["gain_16channel"] == False:
# img /= self.gain
# img.array[img.array > 65535] = 65535
# img.replaceNegative(replace_value=0)
# img.quantize()
# ######################################################################################
# # Output images for calibration pointing
# ######################################################################################
# # Bias output
# if config["ins_effects"]["add_bias"] == True and config["output_setting"]["bias_output"] == True and pointing_type == 'CAL':
# if self.logger is not None:
# self.logger.info(" Output N frame Bias files")
# else:
# print(" Output N frame Bias files", flush=True)
# NBias = int(config["output_setting"]["NBias"])
# for i in range(NBias):
# # BiasCombImg, BiasTag = effects.MakeBiasNcomb(
# # self.npix_x, self.npix_y,
# # bias_level=float(self.bias_level),
# # ncombine=1, read_noise=self.read_noise, gain=1,
# # seed=SeedBiasNonuni+self.chipID,
# # logger=self.logger)
# BiasCombImg = galsim.Image(
# self.npix_x, self.npix_y, init_value=0)
# if config["ins_effects"]["add_bias"] == True:
# biaslevel = self.bias_level
# overscan = biaslevel-2
# elif config["ins_effects"]["add_bias"] == False:
# biaslevel = 0
# overscan = 0
# # Readout noise for Biases is not generated with random seeds. So readout noise for bias images can't be reproduced.
# if config["ins_effects"]["cosmic_ray"] == True:
# if config["ins_effects"]["cray_differ"] == True:
# cr_map, cr_event_num = effects.produceCR_Map(
# xLen=self.npix_x, yLen=self.npix_y,
# exTime=0.01,
# cr_pixelRatio=0.003 *
# (0.01+0.5*self.readout_time)/150.,
# gain=self.gain,
# attachedSizes=self.attachedSizes,
# seed=SeedCosmicRay+pointing_ID*30+self.chipID+1)
# # seed: obj-imaging:+0; bias:+1; dark:+2; flat:+3;
# BiasCombImg += cr_map
# del cr_map
# # Apply Bad lines
# if config["ins_effects"]["add_badcolumns"] == True:
# BiasCombImg = effects.BadColumns(
# BiasCombImg-float(self.bias_level)+5, seed=SeedBadColumns, chipid=self.chipID, logger=self.logger) + float(self.bias_level)-5
# # Non-Linearity for Bias
# if config["ins_effects"]["non_linear"] == True:
# if self.logger is not None:
# self.logger.info(
# " Applying Non-Linearity on the Bias image")
# else:
# print(
# " Applying Non-Linearity on the Bias image", flush=True)
# BiasCombImg = effects.NonLinearity(
# GSImage=BiasCombImg, beta1=5.e-7, beta2=0)
# # START
# pre1 = self.prescan_x # 27
# over1 = self.overscan_x # 71
# pre2 = self.prescan_y # 0 #4
# over2 = self.overscan_y # 84 #80
# # prescan & overscan
# if config["ins_effects"]["add_prescan"] == True:
# chip_utils.log_info(
# msg=" Apply pre/over-scan", logger=self.logger)
# BiasCombImg = chip_utils.AddPreScan(
# GSImage=BiasCombImg, pre1=pre1, pre2=pre2, over1=over1, over2=over2)
# # 1*16 output
# if config["ins_effects"]["format_output"] == True:
# chip_utils.log_info(
# msg=" Apply 1*16 format", logger=self.logger)
# BiasCombImg = chip_utils.formatOutput(GSImage=BiasCombImg)
# self.nsecy = 1
# self.nsecx = 16
# # END
# # Add Bias level
# if config["ins_effects"]["add_bias"] == True:
# if self.logger is not None:
# self.logger.info(
# " Adding Bias level and 16-channel non-uniformity")
# else:
# print(" Adding Bias level and 16-channel non-uniformity")
# BiasCombImg = effects.AddBiasNonUniform16(BiasCombImg,
# bias_level=biaslevel,
# nsecy=self.nsecy, nsecx=self.nsecx,
# seed=SeedBiasNonuni+self.chipID,
# logger=self.logger)
# rng = galsim.UniformDeviate()
# ncombine = 1
# NoiseBias = galsim.GaussianNoise(
# rng=rng, sigma=self.read_noise*ncombine**0.5)
# BiasCombImg.addNoise(NoiseBias)
# BiasCombImg, self.gain_channel = effects.ApplyGainNonUniform16(BiasCombImg, gain=self.gain,
# nsecy=self.nsecy, nsecx=self.nsecx,
# seed=SeedGainNonuni+self.chipID,
# logger=self.logger)
# # BiasCombImg = effects.AddOverscan(
# # BiasCombImg,
# # overscan=float(config["ins_effects"]["bias_level"])-2, gain=self.gain,
# # widthl=27, widthr=27, widtht=8, widthb=8)
# BiasCombImg.replaceNegative(replace_value=0)
# BiasCombImg.quantize()
# BiasCombImg = galsim.ImageUS(BiasCombImg)
# timestamp_obs += 10 * 60
# chip_utils.outputCal(
# chip=self,
# img=BiasCombImg,
# ra_cen=ra_cen,
# dec_cen=dec_cen,
# img_rot=img_rot,
# im_type='BIAS',
# pointing_ID=pointing_ID,
# output_dir=chip_output.subdir,
# exptime=0.0,
# project_cycle=config["project_cycle"],
# run_counter=config["run_counter"],
# timestamp=timestamp_obs)
# del BiasCombImg
# # Export combined (ncombine, Vignetting + PRNU) & single vignetting flat-field file
# if config["ins_effects"]["flat_fielding"] == True and config["output_setting"]["flat_output"] == True and pointing_type == 'CAL':
# if self.logger is not None:
# self.logger.info(" Output N frame Flat-Field files")
# else:
# print(" Output N frame Flat-Field files", flush=True)
# NFlat = int(config["output_setting"]["NFlat"])
# if config["ins_effects"]["add_bias"] == True:
# biaslevel = self.bias_level
# overscan = biaslevel-2
# elif config["ins_effects"]["add_bias"] == False:
# biaslevel = 0
# overscan = 0
# darklevel = self.dark_noise * \
# (self.flat_exptime+0.5*self.readout_time)
# for i in range(NFlat):
# FlatSingle = flat_img * prnu_img + darklevel
# FlatCombImg, FlatTag = effects.MakeFlatNcomb(
# flat_single_image=FlatSingle,
# ncombine=1,
# read_noise=self.read_noise,
# gain=1,
# overscan=overscan,
# biaslevel=0,
# seed_bias=SeedDefective+self.chipID,
# logger=self.logger
# )
# if config["ins_effects"]["cosmic_ray"] == True:
# if config["ins_effects"]["cray_differ"] == True:
# cr_map, cr_event_num = effects.produceCR_Map(
# xLen=self.npix_x, yLen=self.npix_y,
# exTime=self.flat_exptime+0.5*self.readout_time,
# cr_pixelRatio=0.003 *
# (self.flat_exptime+0.5*self.readout_time)/150.,
# gain=self.gain,
# attachedSizes=self.attachedSizes,
# seed=SeedCosmicRay+pointing_ID*30+self.chipID+3)
# # seed: obj-imaging:+0; bias:+1; dark:+2; flat:+3;
# FlatCombImg += cr_map
# del cr_map
# # Add Hot Pixels or/and Dead Pixels
# rgbadpix = Generator(PCG64(int(SeedDefective+self.chipID)))
# badfraction = 5E-5*(rgbadpix.random()*0.5+0.7)
# FlatCombImg = effects.DefectivePixels(
# FlatCombImg, IfHotPix=BoolHotPix, IfDeadPix=BoolDeadPix, fraction=badfraction, seed=SeedDefective+self.chipID, biaslevel=0)
# # Apply Bad lines
# if config["ins_effects"]["add_badcolumns"] == True:
# FlatCombImg = effects.BadColumns(
# FlatCombImg, seed=SeedBadColumns, chipid=self.chipID, logger=self.logger)
# if config["ins_effects"]["non_linear"] == True:
# if self.logger is not None:
# self.logger.info(
# " Applying Non-Linearity on the Flat image")
# else:
# print(
# " Applying Non-Linearity on the Flat image", flush=True)
# FlatCombImg = effects.NonLinearity(
# GSImage=FlatCombImg, beta1=5.e-7, beta2=0)
# # if config["ins_effects"]["cte_trail"] == True:
# # FlatCombImg = effects.CTE_Effect(GSImage=FlatCombImg, threshold=3)
# # START
# pre1 = self.prescan_x # 27
# over1 = self.overscan_x # 71
# pre2 = self.prescan_y # 0 #4
# over2 = self.overscan_y # 84 #80
# if config["ins_effects"]["cte_trail"] == True:
# FlatCombImg = effects.CTE_Effect(GSImage=FlatCombImg, threshold=3)
# START
pre1 = self.prescan_x # 27
over1 = self.overscan_x # 71
pre2 = self.prescan_y # 0 #4
over2 = self.overscan_y # 84 #80
if config["ins_effects"]["cte_trail"] == True:
chip_utils.log_info(
msg=" Apply CTE Effect", logger=self.logger)
# img = effects.CTE_Effect(GSImage=img, threshold=27)
# CTI_modeling
# 2*8 -> 1*16 img-layout
FlatCombImg = chip_utils.formatOutput(GSImage=FlatCombImg)
self.nsecy = 1
self.nsecx = 16
img_arr = FlatCombImg.array
ny, nx = img_arr.shape
dx = int(nx/self.nsecx)
dy = int(ny/self.nsecy)
newimg = galsim.Image(nx, int(ny+over2), init_value=0)
for ichannel in range(16):
print('\n***add CTI effects: pointing-{:} chip-{:} channel-{:}***'.format(
pointing_ID, self.chipID, ichannel+1))
# nx,ny,noverscan,nsp,nmax = 4608,4616,84,3,10
noverscan, nsp, nmax = over2, 3, 10
beta, w, c = 0.478, 84700, 0
t = np.array([0.74, 7.7, 37], dtype=np.float32)
rho_trap = np.array([0.6, 1.6, 1.4], dtype=np.float32)
trap_seeds = np.array(
[0, 1000, 10000], dtype=np.int32) + ichannel + self.chipID*16
release_seed = 50 + ichannel + pointing_ID*30 + self.chipID*16
newimg.array[:, 0+ichannel*dx:dx+ichannel*dx] = CTI_sim(
img_arr[:, 0+ichannel*dx:dx+ichannel*dx], dx, dy, noverscan, nsp, nmax, beta, w, c, t, rho_trap, trap_seeds, release_seed)
newimg.wcs = FlatCombImg.wcs
del FlatCombImg
FlatCombImg = newimg
# 1*16 -> 2*8 img-layout
FlatCombImg = chip_utils.formatRevert(GSImage=FlatCombImg)
self.nsecy = 2
self.nsecx = 8
# prescan & overscan
if config["ins_effects"]["add_prescan"] == True:
chip_utils.log_info(
msg=" Apply pre/over-scan", logger=self.logger)
if config["ins_effects"]["cte_trail"] == False:
FlatCombImg = chip_utils.AddPreScan(
GSImage=FlatCombImg, pre1=pre1, pre2=pre2, over1=over1, over2=over2)
if config["ins_effects"]["cte_trail"] == True:
FlatCombImg = chip_utils.AddPreScan(
GSImage=FlatCombImg, pre1=pre1, pre2=pre2, over1=over1, over2=0)
# 1*16 output
if config["ins_effects"]["format_output"] == True:
chip_utils.log_info(
msg=" Apply 1*16 format", logger=self.logger)
FlatCombImg = chip_utils.formatOutput(GSImage=FlatCombImg)
self.nsecy = 1
self.nsecx = 16
# END
# Add Bias level
if config["ins_effects"]["add_bias"] == True:
if self.logger is not None:
self.logger.info(
" Adding Bias level and 16-channel non-uniformity")
else:
print(" Adding Bias level and 16-channel non-uniformity")
# img += float(config["ins_effects"]["bias_level"])
FlatCombImg = effects.AddBiasNonUniform16(FlatCombImg,
bias_level=biaslevel,
nsecy=self.nsecy, nsecx=self.nsecx,
seed=SeedBiasNonuni+self.chipID,
logger=self.logger)
# Add Read-out Noise
if config["ins_effects"]["add_readout"] == True:
seed = int(config["random_seeds"]["seed_readout"]
) + pointing_ID*30 + self.chipID + 3
rng_readout = galsim.BaseDeviate(seed)
readout_noise = galsim.GaussianNoise(
rng=rng_readout, sigma=self.read_noise)
FlatCombImg.addNoise(readout_noise)
FlatCombImg, self.gain_channel = effects.ApplyGainNonUniform16(FlatCombImg, gain=self.gain,
nsecy=self.nsecy, nsecx=self.nsecx,
seed=SeedGainNonuni+self.chipID,
logger=self.logger)
# FlatCombImg = effects.AddOverscan(FlatCombImg, overscan=overscan, gain=self.gain, widthl=27, widthr=27, widtht=8, widthb=8)
FlatCombImg.replaceNegative(replace_value=0)
FlatCombImg.quantize()
FlatCombImg = galsim.ImageUS(FlatCombImg)
timestamp_obs += 10 * 60
chip_utils.outputCal(
chip=self,
img=FlatCombImg,
ra_cen=ra_cen,
dec_cen=dec_cen,
img_rot=img_rot,
im_type='FLAT',
pointing_ID=pointing_ID,
output_dir=chip_output.subdir,
exptime=self.flat_exptime,
project_cycle=config["project_cycle"],
run_counter=config["run_counter"],
timestamp=timestamp_obs)
del FlatCombImg, FlatSingle, prnu_img
# flat_img.replaceNegative(replace_value=0)
# flat_img.quantize()
# galsim.ImageUS(flat_img).write("%s/FlatImg_Vignette_%s.fits" % (chip_output.subdir, self.chipID))
del flat_img
# Export Dark current images
if config["ins_effects"]["add_dark"] == True and config["output_setting"]["dark_output"] == True and pointing_type == 'CAL':
if self.logger is not None:
self.logger.info(" Output N frame Dark Current files")
else:
print(" Output N frame Dark Current files", flush=True)
NDark = int(config["output_setting"]["NDark"])
if config["ins_effects"]["add_bias"] == True:
biaslevel = self.bias_level
overscan = biaslevel-2
elif config["ins_effects"]["add_bias"] == False:
biaslevel = 0
overscan = 0
for i in range(NDark):
DarkCombImg, DarkTag = effects.MakeDarkNcomb(
self.npix_x, self.npix_y,
overscan=overscan, bias_level=0, darkpsec=0.02, exptime=self.dark_exptime+0.5*self.readout_time,
ncombine=1, read_noise=self.read_noise,
gain=1, seed_bias=SeedBiasNonuni+self.chipID,
logger=self.logger)
if config["ins_effects"]["cosmic_ray"] == True:
if config["ins_effects"]["cray_differ"] == True:
cr_map, cr_event_num = effects.produceCR_Map(
xLen=self.npix_x, yLen=self.npix_y,
exTime=self.dark_exptime+0.5*self.readout_time,
cr_pixelRatio=0.003 *
(self.dark_exptime+0.5*self.readout_time)/150.,
gain=self.gain,
attachedSizes=self.attachedSizes,
seed=SeedCosmicRay+pointing_ID*30+self.chipID+2)
# seed: obj-imaging:+0; bias:+1; dark:+2; flat:+3;
DarkCombImg += cr_map
cr_map[cr_map > 65535] = 65535
cr_map[cr_map < 0] = 0
crmap_gsimg = galsim.Image(cr_map, dtype=np.uint16)
del cr_map
# START
# prescan & overscan
if config["ins_effects"]["add_prescan"] == True:
chip_utils.log_info(
msg=" Apply pre/over-scan", logger=self.logger)
crmap_gsimg = chip_utils.AddPreScan(
GSImage=crmap_gsimg, pre1=pre1, pre2=pre2, over1=over1, over2=over2)
# 1*16 output
if config["ins_effects"]["format_output"] == True:
chip_utils.log_info(
msg=" Apply 1*16 format", logger=self.logger)
crmap_gsimg = chip_utils.formatOutput(
GSImage=crmap_gsimg)
self.nsecy = 1
self.nsecx = 16
# END
chip_utils.outputCal(
chip=self,
img=crmap_gsimg,
ra_cen=ra_cen,
dec_cen=dec_cen,
img_rot=img_rot,
im_type='CRD',
pointing_ID=pointing_ID,
output_dir=chip_output.subdir,
exptime=self.dark_exptime,
project_cycle=config["project_cycle"],
run_counter=config["run_counter"],
timestamp=timestamp_obs)
del crmap_gsimg
# Add Hot Pixels or/and Dead Pixels
rgbadpix = Generator(PCG64(int(SeedDefective+self.chipID)))
badfraction = 5E-5*(rgbadpix.random()*0.5+0.7)
DarkCombImg = effects.DefectivePixels(
DarkCombImg, IfHotPix=BoolHotPix, IfDeadPix=BoolDeadPix, fraction=badfraction, seed=SeedDefective+self.chipID, biaslevel=0)
# Apply Bad lines
if config["ins_effects"]["add_badcolumns"] == True:
DarkCombImg = effects.BadColumns(
DarkCombImg, seed=SeedBadColumns, chipid=self.chipID, logger=self.logger)
# Non-Linearity for Dark
if config["ins_effects"]["non_linear"] == True:
if self.logger is not None:
self.logger.info(
" Applying Non-Linearity on the Dark image")
else:
print(
" Applying Non-Linearity on the Dark image", flush=True)
DarkCombImg = effects.NonLinearity(
GSImage=DarkCombImg, beta1=5.e-7, beta2=0)
# chip_utils.log_info(
# msg=" Apply CTE Effect", logger=self.logger)
# # img = effects.CTE_Effect(GSImage=img, threshold=27)
# # CTI_modeling
# # 2*8 -> 1*16 img-layout
# FlatCombImg = chip_utils.formatOutput(GSImage=FlatCombImg)
# self.nsecy = 1
# self.nsecx = 16
# img_arr = FlatCombImg.array
# ny, nx = img_arr.shape
# dx = int(nx/self.nsecx)
# dy = int(ny/self.nsecy)
# newimg = galsim.Image(nx, int(ny+over2), init_value=0)
# for ichannel in range(16):
# print('\n***add CTI effects: pointing-{:} chip-{:} channel-{:}***'.format(
# pointing_ID, self.chipID, ichannel+1))
# # nx,ny,noverscan,nsp,nmax = 4608,4616,84,3,10
# noverscan, nsp, nmax = over2, 3, 10
# beta, w, c = 0.478, 84700, 0
# t = np.array([0.74, 7.7, 37], dtype=np.float32)
# rho_trap = np.array([0.6, 1.6, 1.4], dtype=np.float32)
# trap_seeds = np.array(
# [0, 1000, 10000], dtype=np.int32) + ichannel + self.chipID*16
# release_seed = 50 + ichannel + pointing_ID*30 + self.chipID*16
# newimg.array[:, 0+ichannel*dx:dx+ichannel*dx] = CTI_sim(
# img_arr[:, 0+ichannel*dx:dx+ichannel*dx], dx, dy, noverscan, nsp, nmax, beta, w, c, t, rho_trap, trap_seeds, release_seed)
# newimg.wcs = FlatCombImg.wcs
# del FlatCombImg
# FlatCombImg = newimg
# # 1*16 -> 2*8 img-layout
# FlatCombImg = chip_utils.formatRevert(GSImage=FlatCombImg)
# self.nsecy = 2
# self.nsecx = 8
# # prescan & overscan
# if config["ins_effects"]["add_prescan"] == True:
# chip_utils.log_info(
# msg=" Apply pre/over-scan", logger=self.logger)
# if config["ins_effects"]["cte_trail"] == False:
# FlatCombImg = chip_utils.AddPreScan(
# GSImage=FlatCombImg, pre1=pre1, pre2=pre2, over1=over1, over2=over2)
# if config["ins_effects"]["cte_trail"] == True:
# FlatCombImg = chip_utils.AddPreScan(
# GSImage=FlatCombImg, pre1=pre1, pre2=pre2, over1=over1, over2=0)
# # 1*16 output
# if config["ins_effects"]["format_output"] == True:
# chip_utils.log_info(
# msg=" Apply 1*16 format", logger=self.logger)
# FlatCombImg = chip_utils.formatOutput(GSImage=FlatCombImg)
# self.nsecy = 1
# self.nsecx = 16
# # END
# # Add Bias level
# if config["ins_effects"]["add_bias"] == True:
# if self.logger is not None:
# self.logger.info(
# " Adding Bias level and 16-channel non-uniformity")
# else:
# print(" Adding Bias level and 16-channel non-uniformity")
# # img += float(config["ins_effects"]["bias_level"])
# FlatCombImg = effects.AddBiasNonUniform16(FlatCombImg,
# bias_level=biaslevel,
# nsecy=self.nsecy, nsecx=self.nsecx,
# seed=SeedBiasNonuni+self.chipID,
# logger=self.logger)
# # Add Read-out Noise
# if config["ins_effects"]["add_readout"] == True:
# seed = int(config["random_seeds"]["seed_readout"]
# ) + pointing_ID*30 + self.chipID + 3
# rng_readout = galsim.BaseDeviate(seed)
# readout_noise = galsim.GaussianNoise(
# rng=rng_readout, sigma=self.read_noise)
# FlatCombImg.addNoise(readout_noise)
# FlatCombImg, self.gain_channel = effects.ApplyGainNonUniform16(FlatCombImg, gain=self.gain,
# nsecy=self.nsecy, nsecx=self.nsecx,
# seed=SeedGainNonuni+self.chipID,
# logger=self.logger)
# # FlatCombImg = effects.AddOverscan(FlatCombImg, overscan=overscan, gain=self.gain, widthl=27, widthr=27, widtht=8, widthb=8)
# FlatCombImg.replaceNegative(replace_value=0)
# FlatCombImg.quantize()
# FlatCombImg = galsim.ImageUS(FlatCombImg)
# timestamp_obs += 10 * 60
# chip_utils.outputCal(
# chip=self,
# img=FlatCombImg,
# ra_cen=ra_cen,
# dec_cen=dec_cen,
# img_rot=img_rot,
# im_type='FLAT',
# pointing_ID=pointing_ID,
# output_dir=chip_output.subdir,
# exptime=self.flat_exptime,
# project_cycle=config["project_cycle"],
# run_counter=config["run_counter"],
# timestamp=timestamp_obs)
# del FlatCombImg, FlatSingle, prnu_img
# # flat_img.replaceNegative(replace_value=0)
# # flat_img.quantize()
# # galsim.ImageUS(flat_img).write("%s/FlatImg_Vignette_%s.fits" % (chip_output.subdir, self.chipID))
# del flat_img
# # Export Dark current images
# if config["ins_effects"]["add_dark"] == True and config["output_setting"]["dark_output"] == True and pointing_type == 'CAL':
# if self.logger is not None:
# self.logger.info(" Output N frame Dark Current files")
# else:
# print(" Output N frame Dark Current files", flush=True)
# NDark = int(config["output_setting"]["NDark"])
# if config["ins_effects"]["add_bias"] == True:
# biaslevel = self.bias_level
# overscan = biaslevel-2
# elif config["ins_effects"]["add_bias"] == False:
# biaslevel = 0
# overscan = 0
# for i in range(NDark):
# DarkCombImg, DarkTag = effects.MakeDarkNcomb(
# self.npix_x, self.npix_y,
# overscan=overscan, bias_level=0, darkpsec=0.02, exptime=self.dark_exptime+0.5*self.readout_time,
# ncombine=1, read_noise=self.read_noise,
# gain=1, seed_bias=SeedBiasNonuni+self.chipID,
# logger=self.logger)
# if config["ins_effects"]["cosmic_ray"] == True:
# if config["ins_effects"]["cray_differ"] == True:
# cr_map, cr_event_num = effects.produceCR_Map(
# xLen=self.npix_x, yLen=self.npix_y,
# exTime=self.dark_exptime+0.5*self.readout_time,
# cr_pixelRatio=0.003 *
# (self.dark_exptime+0.5*self.readout_time)/150.,
# gain=self.gain,
# attachedSizes=self.attachedSizes,
# seed=SeedCosmicRay+pointing_ID*30+self.chipID+2)
# # seed: obj-imaging:+0; bias:+1; dark:+2; flat:+3;
# DarkCombImg += cr_map
# cr_map[cr_map > 65535] = 65535
# cr_map[cr_map < 0] = 0
# crmap_gsimg = galsim.Image(cr_map, dtype=np.uint16)
# del cr_map
# # START
# # prescan & overscan
# if config["ins_effects"]["add_prescan"] == True:
# chip_utils.log_info(
# msg=" Apply pre/over-scan", logger=self.logger)
# crmap_gsimg = chip_utils.AddPreScan(
# GSImage=crmap_gsimg, pre1=pre1, pre2=pre2, over1=over1, over2=over2)
# # 1*16 output
# if config["ins_effects"]["format_output"] == True:
# chip_utils.log_info(
# msg=" Apply 1*16 format", logger=self.logger)
# crmap_gsimg = chip_utils.formatOutput(
# GSImage=crmap_gsimg)
# self.nsecy = 1
# self.nsecx = 16
# # END
# chip_utils.outputCal(
# chip=self,
# img=crmap_gsimg,
# ra_cen=ra_cen,
# dec_cen=dec_cen,
# img_rot=img_rot,
# im_type='CRD',
# pointing_ID=pointing_ID,
# output_dir=chip_output.subdir,
# exptime=self.dark_exptime,
# project_cycle=config["project_cycle"],
# run_counter=config["run_counter"],
# timestamp=timestamp_obs)
# del crmap_gsimg
# # Add Hot Pixels or/and Dead Pixels
# rgbadpix = Generator(PCG64(int(SeedDefective+self.chipID)))
# badfraction = 5E-5*(rgbadpix.random()*0.5+0.7)
# DarkCombImg = effects.DefectivePixels(
# DarkCombImg, IfHotPix=BoolHotPix, IfDeadPix=BoolDeadPix, fraction=badfraction, seed=SeedDefective+self.chipID, biaslevel=0)
# # Apply Bad lines
# if config["ins_effects"]["add_badcolumns"] == True:
# DarkCombImg = effects.BadColumns(
# DarkCombImg, seed=SeedBadColumns, chipid=self.chipID, logger=self.logger)
# # Non-Linearity for Dark
# if config["ins_effects"]["non_linear"] == True:
# if self.logger is not None:
# self.logger.info(
# " Applying Non-Linearity on the Dark image")
# else:
# print(
# " Applying Non-Linearity on the Dark image", flush=True)
# DarkCombImg = effects.NonLinearity(
# GSImage=DarkCombImg, beta1=5.e-7, beta2=0)
# # if config["ins_effects"]["cte_trail"] == True:
# # DarkCombImg = effects.CTE_Effect(GSImage=DarkCombImg, threshold=3)
# # START
# pre1 = self.prescan_x # 27
# over1 = self.overscan_x # 71
# pre2 = self.prescan_y # 0 #4
# over2 = self.overscan_y # 84 #80
# if config["ins_effects"]["cte_trail"] == True:
# DarkCombImg = effects.CTE_Effect(GSImage=DarkCombImg, threshold=3)
# START
pre1 = self.prescan_x # 27
over1 = self.overscan_x # 71
pre2 = self.prescan_y # 0 #4
over2 = self.overscan_y # 84 #80
if config["ins_effects"]["cte_trail"] == True:
chip_utils.log_info(
msg=" Apply CTE Effect", logger=self.logger)
# img = effects.CTE_Effect(GSImage=img, threshold=27)
# CTI_modeling
# 2*8 -> 1*16 img-layout
DarkCombImg = chip_utils.formatOutput(GSImage=DarkCombImg)
self.nsecy = 1
self.nsecx = 16
img_arr = DarkCombImg.array
ny, nx = img_arr.shape
dx = int(nx/self.nsecx)
dy = int(ny/self.nsecy)
newimg = galsim.Image(nx, int(ny+over2), init_value=0)
for ichannel in range(16):
print('\n***add CTI effects: pointing-{:} chip-{:} channel-{:}***'.format(
pointing_ID, self.chipID, ichannel+1))
# nx,ny,noverscan,nsp,nmax = 4608,4616,84,3,10
noverscan, nsp, nmax = over2, 3, 10
beta, w, c = 0.478, 84700, 0
t = np.array([0.74, 7.7, 37], dtype=np.float32)
rho_trap = np.array([0.6, 1.6, 1.4], dtype=np.float32)
trap_seeds = np.array(
[0, 1000, 10000], dtype=np.int32) + ichannel + self.chipID*16
release_seed = 50 + ichannel + pointing_ID*30 + self.chipID*16
newimg.array[:, 0+ichannel*dx:dx+ichannel*dx] = CTI_sim(
img_arr[:, 0+ichannel*dx:dx+ichannel*dx], dx, dy, noverscan, nsp, nmax, beta, w, c, t, rho_trap, trap_seeds, release_seed)
newimg.wcs = DarkCombImg.wcs
del DarkCombImg
DarkCombImg = newimg
# 1*16 -> 2*8 img-layout
DarkCombImg = chip_utils.formatRevert(GSImage=DarkCombImg)
self.nsecy = 2
self.nsecx = 8
# chip_utils.log_info(
# msg=" Apply CTE Effect", logger=self.logger)
# # img = effects.CTE_Effect(GSImage=img, threshold=27)
# # CTI_modeling
# # 2*8 -> 1*16 img-layout
# DarkCombImg = chip_utils.formatOutput(GSImage=DarkCombImg)
# self.nsecy = 1
# self.nsecx = 16
# img_arr = DarkCombImg.array
# ny, nx = img_arr.shape
# dx = int(nx/self.nsecx)
# dy = int(ny/self.nsecy)
# newimg = galsim.Image(nx, int(ny+over2), init_value=0)
# for ichannel in range(16):
# print('\n***add CTI effects: pointing-{:} chip-{:} channel-{:}***'.format(
# pointing_ID, self.chipID, ichannel+1))
# # nx,ny,noverscan,nsp,nmax = 4608,4616,84,3,10
# noverscan, nsp, nmax = over2, 3, 10
# beta, w, c = 0.478, 84700, 0
# t = np.array([0.74, 7.7, 37], dtype=np.float32)
# rho_trap = np.array([0.6, 1.6, 1.4], dtype=np.float32)
# trap_seeds = np.array(
# [0, 1000, 10000], dtype=np.int32) + ichannel + self.chipID*16
# release_seed = 50 + ichannel + pointing_ID*30 + self.chipID*16
# newimg.array[:, 0+ichannel*dx:dx+ichannel*dx] = CTI_sim(
# img_arr[:, 0+ichannel*dx:dx+ichannel*dx], dx, dy, noverscan, nsp, nmax, beta, w, c, t, rho_trap, trap_seeds, release_seed)
# newimg.wcs = DarkCombImg.wcs
# del DarkCombImg
# DarkCombImg = newimg
# # 1*16 -> 2*8 img-layout
# DarkCombImg = chip_utils.formatRevert(GSImage=DarkCombImg)
# self.nsecy = 2
# self.nsecx = 8
# # prescan & overscan
# if config["ins_effects"]["add_prescan"] == True:
# chip_utils.log_info(
# msg=" Apply pre/over-scan", logger=self.logger)
# if config["ins_effects"]["cte_trail"] == False:
# DarkCombImg = chip_utils.AddPreScan(
# GSImage=DarkCombImg, pre1=pre1, pre2=pre2, over1=over1, over2=over2)
# if config["ins_effects"]["cte_trail"] == True:
# DarkCombImg = chip_utils.AddPreScan(
# GSImage=DarkCombImg, pre1=pre1, pre2=pre2, over1=over1, over2=0)
# # 1*16 output
# if config["ins_effects"]["format_output"] == True:
# chip_utils.log_info(
# msg=" Apply 1*16 format", logger=self.logger)
# DarkCombImg = chip_utils.formatOutput(GSImage=DarkCombImg)
# self.nsecy = 1
# self.nsecx = 16
# # END
# # Add Bias level
# if config["ins_effects"]["add_bias"] == True:
# if self.logger is not None:
# self.logger.info(
# " Adding Bias level and 16-channel non-uniformity")
# else:
# print(" Adding Bias level and 16-channel non-uniformity")
# # img += float(config["ins_effects"]["bias_level"])
# DarkCombImg = effects.AddBiasNonUniform16(DarkCombImg,
# bias_level=biaslevel,
# nsecy=self.nsecy, nsecx=self.nsecx,
# seed=SeedBiasNonuni+self.chipID,
# logger=self.logger)
# prescan & overscan
if config["ins_effects"]["add_prescan"] == True:
chip_utils.log_info(
msg=" Apply pre/over-scan", logger=self.logger)
if config["ins_effects"]["cte_trail"] == False:
DarkCombImg = chip_utils.AddPreScan(
GSImage=DarkCombImg, pre1=pre1, pre2=pre2, over1=over1, over2=over2)
if config["ins_effects"]["cte_trail"] == True:
DarkCombImg = chip_utils.AddPreScan(
GSImage=DarkCombImg, pre1=pre1, pre2=pre2, over1=over1, over2=0)
# 1*16 output
if config["ins_effects"]["format_output"] == True:
chip_utils.log_info(
msg=" Apply 1*16 format", logger=self.logger)
DarkCombImg = chip_utils.formatOutput(GSImage=DarkCombImg)
self.nsecy = 1
self.nsecx = 16
# END
# Add Bias level
if config["ins_effects"]["add_bias"] == True:
if self.logger is not None:
self.logger.info(
" Adding Bias level and 16-channel non-uniformity")
else:
print(" Adding Bias level and 16-channel non-uniformity")
# img += float(config["ins_effects"]["bias_level"])
DarkCombImg = effects.AddBiasNonUniform16(DarkCombImg,
bias_level=biaslevel,
nsecy=self.nsecy, nsecx=self.nsecx,
seed=SeedBiasNonuni+self.chipID,
logger=self.logger)
# Add Read-out Noise
if config["ins_effects"]["add_readout"] == True:
seed = int(config["random_seeds"]["seed_readout"]
) + pointing_ID*30 + self.chipID + 2
rng_readout = galsim.BaseDeviate(seed)
readout_noise = galsim.GaussianNoise(
rng=rng_readout, sigma=self.read_noise)
DarkCombImg.addNoise(readout_noise)
DarkCombImg, self.gain_channel = effects.ApplyGainNonUniform16(
DarkCombImg, gain=self.gain,
nsecy=self.nsecy, nsecx=self.nsecx,
seed=SeedGainNonuni+self.chipID,
logger=self.logger)
# DarkCombImg = effects.AddOverscan(
# DarkCombImg,
# overscan=overscan, gain=self.gain,
# widthl=27, widthr=27, widtht=8, widthb=8)
DarkCombImg.replaceNegative(replace_value=0)
DarkCombImg.quantize()
DarkCombImg = galsim.ImageUS(DarkCombImg)
timestamp_obs += 10 * 60
chip_utils.outputCal(
chip=self,
img=DarkCombImg,
ra_cen=ra_cen,
dec_cen=dec_cen,
img_rot=img_rot,
im_type='DARK',
pointing_ID=pointing_ID,
output_dir=chip_output.subdir,
exptime=self.dark_exptime,
project_cycle=config["project_cycle"],
run_counter=config["run_counter"],
timestamp=timestamp_obs)
del DarkCombImg
# img = galsim.ImageUS(img)
# # 16 output channel, with each a single image file
# if config["ins_effects"]["readout16"] == True:
# print(" 16 Output Channel simulation")
# for coli in [0, 1]:
# for rowi in range(8):
# sub_img = effects.readout16(
# GSImage=img,
# rowi=rowi,
# coli=coli,
# overscan_value=self.overscan)
# rowcoltag = str(rowi) + str(coli)
# img_name_root = chip_output.img_name.split(".")[0]
# sub_img.write("%s/%s_%s.fits" % (chip_output.subdir, img_name_root, rowcoltag))
# del sub_img
return img
# # Add Read-out Noise
# if config["ins_effects"]["add_readout"] == True:
# seed = int(config["random_seeds"]["seed_readout"]
# ) + pointing_ID*30 + self.chipID + 2
# rng_readout = galsim.BaseDeviate(seed)
# readout_noise = galsim.GaussianNoise(
# rng=rng_readout, sigma=self.read_noise)
# DarkCombImg.addNoise(readout_noise)
# DarkCombImg, self.gain_channel = effects.ApplyGainNonUniform16(
# DarkCombImg, gain=self.gain,
# nsecy=self.nsecy, nsecx=self.nsecx,
# seed=SeedGainNonuni+self.chipID,
# logger=self.logger)
# # DarkCombImg = effects.AddOverscan(
# # DarkCombImg,
# # overscan=overscan, gain=self.gain,
# # widthl=27, widthr=27, widtht=8, widthb=8)
# DarkCombImg.replaceNegative(replace_value=0)
# DarkCombImg.quantize()
# DarkCombImg = galsim.ImageUS(DarkCombImg)
# timestamp_obs += 10 * 60
# chip_utils.outputCal(
# chip=self,
# img=DarkCombImg,
# ra_cen=ra_cen,
# dec_cen=dec_cen,
# img_rot=img_rot,
# im_type='DARK',
# pointing_ID=pointing_ID,
# output_dir=chip_output.subdir,
# exptime=self.dark_exptime,
# project_cycle=config["project_cycle"],
# run_counter=config["run_counter"],
# timestamp=timestamp_obs)
# del DarkCombImg
# # img = galsim.ImageUS(img)
# # # 16 output channel, with each a single image file
# # if config["ins_effects"]["readout16"] == True:
# # print(" 16 Output Channel simulation")
# # for coli in [0, 1]:
# # for rowi in range(8):
# # sub_img = effects.readout16(
# # GSImage=img,
# # rowi=rowi,
# # coli=coli,
# # overscan_value=self.overscan)
# # rowcoltag = str(rowi) + str(coli)
# # img_name_root = chip_output.img_name.split(".")[0]
# # sub_img.write("%s/%s_%s.fits" % (chip_output.subdir, img_name_root, rowcoltag))
# # del sub_img
# return img
ObservationSim/Instrument/FocalPlane.py
View file @ fe70f366
import galsim
import numpy as np
class FocalPlane(object):
def __init__(self, config=None, chip_list=None, survey_type='Photometric', bad_chips=None):
def __init__(self, chip_list=None, survey_type='Photometric', bad_chips=None):
"""Get the focal plane layout
"""
self.nchips = 42
......@@ -39,18 +40,18 @@ class FocalPlane(object):
for i in range(1, 31):
self.ignore_chips.append(i)
if config is not None:
self.nchip_x = config["nchip_x"]
self.nchip_y = config["nchip_y"]
self.npix_tot_x = config["npix_tot_x"]
self.npix_tot_y = config["npix_tot_y"]
self.npix_gap_x = config["npix_gap_x"]
self.npix_gap_y = config["npix_gap_y"]
if "chipLabelIDs" in config:
self.chipLabelIDs = config["chipLabelIDs"]
if "bad_chips" in config:
self.bad_chips = config["bad_chips"]
else:
# if config is not None:
# self.nchip_x = config["nchip_x"]
# self.nchip_y = config["nchip_y"]
# self.npix_tot_x = config["npix_tot_x"]
# self.npix_tot_y = config["npix_tot_y"]
# self.npix_gap_x = config["npix_gap_x"]
# self.npix_gap_y = config["npix_gap_y"]
# if "chipLabelIDs" in config:
# self.chipLabelIDs = config["chipLabelIDs"]
# if "bad_chips" in config:
# self.bad_chips = config["bad_chips"]
# else:
self.nchip_x = 6
self.nchip_y = 5
self.npix_tot_x = 59516
......@@ -65,7 +66,7 @@ class FocalPlane(object):
self.cen_pix_y = 0
def getChipLabel(self, chipID):
return str("0%d"%chipID)[-2:]
return str("0%d" % chipID)[-2:]
def isBadChip(self, chipID):
"""Check if chip #(chipID) on the focal plane is bad or not
......@@ -89,7 +90,8 @@ class FocalPlane(object):
WCS of the focal plane
"""
if logger is not None:
logger.info(" Construct the wcs of the entire image mosaic using Gnomonic/TAN projection")
logger.info(
" Construct the wcs of the entire image mosaic using Gnomonic/TAN projection")
if (xcen == None) or (ycen == None):
xcen = self.cen_pix_x
ycen = self.cen_pix_y
......@@ -104,7 +106,8 @@ class FocalPlane(object):
dvdy = -np.cos(img_rot.rad) * pix_scale
moscen = galsim.PositionD(x=xcen, y=ycen)
sky_center = galsim.CelestialCoord(ra=ra*galsim.degrees, dec=dec*galsim.degrees)
sky_center = galsim.CelestialCoord(
ra=ra*galsim.degrees, dec=dec*galsim.degrees)
affine = galsim.AffineTransform(dudx, dudy, dvdx, dvdy, origin=moscen)
WCS = galsim.TanWCS(affine, sky_center, units=galsim.arcsec)
......@@ -115,10 +118,10 @@ class FocalPlane(object):
The sky coverage of an area
"""
r2d = 180.0/np.pi
s1 = wcs.toWorld(galsim.PositionD(x0,y0))
s2 = wcs.toWorld(galsim.PositionD(x0,y1))
s3 = wcs.toWorld(galsim.PositionD(x1,y0))
s4 = wcs.toWorld(galsim.PositionD(x1,y1))
s1 = wcs.toWorld(galsim.PositionD(x0, y0))
s2 = wcs.toWorld(galsim.PositionD(x0, y1))
s3 = wcs.toWorld(galsim.PositionD(x1, y0))
s4 = wcs.toWorld(galsim.PositionD(x1, y1))
ra = [s1.ra.rad*r2d, s2.ra.rad*r2d, s3.ra.rad*r2d, s4.ra.rad*r2d]
dec = [s1.dec.rad*r2d, s2.dec.rad*r2d, s3.dec.rad*r2d, s4.dec.rad*r2d]
......
tests/test_focalplane.py 0 → 100644
View file @ fe70f366
import unittest
import os
import galsim
from ObservationSim.Instrument import FocalPlane, Chip
class TestFocalPlane(unittest.TestCase):
def __init__(self, methodName='runTest'):
super(TestFocalPlane, self).__init__(methodName)
self.dataPath = os.path.join(
os.getenv('UNIT_TEST_DATA_ROOT'), 'csst_msc_sim/csst_fz_msc')
self.focal_plane = FocalPlane(
chip_list=['8'])
self.assertTrue(self.focal_plane.cen_pix_x == 0)
self.assertTrue(self.focal_plane.cen_pix_y == 0)
test_focal_plane_phot = FocalPlane(survey_type='Photometric')
test_focal_plane_spec = FocalPlane(survey_type='Spectroscopic')
test_focal_plane_FGS = FocalPlane(survey_type='FGS')
test_focal_plane_bad_chips = FocalPlane(bad_chips=["1"])
def test_fp_method(self):
wcs = self.focal_plane.getTanWCS(
192.8595, 0., 0.*galsim.degrees, 0.0074)
sky_coverage = self.focal_plane.getSkyCoverage(
wcs, x0=-1, x1=0, y0=-1, y1=0)
print(sky_coverage.area())
self.assertTrue(abs(sky_coverage.area() - 0.0074**2/(3600.**2)) < 1e13)
if __name__ == '__main_':
unittest.main()
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment