add the simulaiton of flat-field cube for slitless spectroscopic sim (98d94696) · Commits · Liu Dezi / csst_msc_sim

ObservationSim/Instrument/Chip/Chip.py

+19 −0

Original line number	Diff line number	Diff line
		@@ -82,6 +82,8 @@ class Chip(FocalPlane):
		else:
		self.sensor = galsim.Sensor()

		self.flat_cube = None # for spectroscopic flat field cube simulation

		# def _getChipRowCol(self):
		# self.rowID = (self.chipID - 1) // 5 + 1
		# self.colID = (self.chipID - 1) % 5 + 1
		@@ -846,3 +848,20 @@ class Chip(FocalPlane):
		# sub_img.write("%s/%s_%s.fits" % (chip_output.subdir, img_name_root, rowcoltag))
		# del sub_img
		return img

		def loadSLSFLATCUBE(self, flat_fn='flat_cube.fits'):
		from astropy.io import fits
		try:
		with pkg_resources.files('ObservationSim.Instrument.data').joinpath(flat_fn) as data_path:
		flat_fits = fits.open(data_path, ignore_missing_simple=True)
		except AttributeError:
		with pkg_resources.path('ObservationSim.Instrument.data', flat_fn) as data_path:
		flat_fits = fits.open(data_path, ignore_missing_simple=True)

		fl = len(flat_fits)
		fl_sh = flat_fits[0].data.shape
		assert fl == 4, 'FLAT Field Cube is Not 4 layess!!!!!!!'
		self.flat_cube = np.zeros([fl, fl_sh[0], fl_sh[1]])
		for i in np.arange(0, fl, 1):
		self.flat_cube[i] = flat_fits[i].data

ObservationSim/Instrument/data/flatCube/flat_cube.fits

0 → 100644

+2.54 GiB

File added.

No diff preview for this file type.

View file

ObservationSim/MockObject/Galaxy.py

+11 −4

Original line number	Diff line number	Diff line
		@@ -271,6 +271,9 @@ class Galaxy(MockObject):
		folding_threshold = 5.e-3
		gsp = galsim.GSParams(folding_threshold=folding_threshold)
		# nphotons_sum = 0

		flat_cube = chip.flat_cube

		xOrderSigPlus = {'A':1.3909419820029296,'B':1.4760376591236062,'C':4.035447379743442,'D':5.5684364343742825,'E':16.260021029735388}
		grating_split_pos_chip = 0 + grating_split_pos
		for i in range(len(bandpass_list)):
		@@ -319,7 +322,8 @@ class Galaxy(MockObject):
		tar_spec=normalSED,
		band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10,
		conf=chip.sls_conf[0],
		isAlongY=0)
		isAlongY=0,
		flat_cube=flat_cube)

		self.addSLStoChipImage(sdp=sdp_p1, chip=chip, xOrderSigPlus = xOrderSigPlus, local_wcs=real_wcs_local)

		@@ -334,7 +338,8 @@ class Galaxy(MockObject):
		tar_spec=normalSED,
		band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10,
		conf=chip.sls_conf[1],
		isAlongY=0)
		isAlongY=0,
		flat_cube=flat_cube)

		self.addSLStoChipImage(sdp=sdp_p2, chip=chip, xOrderSigPlus = xOrderSigPlus, local_wcs=real_wcs_local)

		@@ -346,7 +351,8 @@ class Galaxy(MockObject):
		tar_spec=normalSED,
		band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10,
		conf=chip.sls_conf[1],
		isAlongY=0)
		isAlongY=0,
		flat_cube=flat_cube)
		self.addSLStoChipImage(sdp=sdp, chip=chip, xOrderSigPlus = xOrderSigPlus, local_wcs=real_wcs_local)
		del sdp
		elif grating_split_pos_chip>=gal_end[1]:
		@@ -355,7 +361,8 @@ class Galaxy(MockObject):
		tar_spec=normalSED,
		band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10,
		conf=chip.sls_conf[0],
		isAlongY=0)
		isAlongY=0,
		flat_cube=flat_cube)
		self.addSLStoChipImage(sdp=sdp, chip=chip, xOrderSigPlus = xOrderSigPlus, local_wcs=real_wcs_local)
		del sdp

ObservationSim/MockObject/MockObject.py

+53 −44

Original line number	Diff line number	Diff line
		@@ -4,9 +4,11 @@ import astropy.constants as cons
		from astropy.table import Table

		from ObservationSim.MockObject._util import magToFlux, VC_A, convolveGaussXorders
		from ObservationSim.MockObject._util import integrate_sed_bandpass, getNormFactorForSpecWithABMAG, getObservedSED, getABMAG
		from ObservationSim.MockObject._util import integrate_sed_bandpass, getNormFactorForSpecWithABMAG, getObservedSED, \
		getABMAG
		from ObservationSim.MockObject.SpecDisperser import SpecDisperser


		class MockObject(object):
		def __init__(self, param, logger=None):
		self.param = param
		@@ -110,7 +112,6 @@ class MockObject(object):

		return realPos


		def drawObject(self, img, final, flux=None, filt=None, tel=None, exptime=150.):
		""" Draw (point like) object on img.
		Should be overided for extended source, e.g. galaxy...
		@@ -139,7 +140,8 @@ class MockObject(object):
		img[bounds] += stamp[bounds]
		return img, stamp, isUpdated

		def drawObj_multiband(self, tel, pos_img, psf_model, bandpass_list, filt, chip, nphotons_tot=None, g1=0, g2=0, exptime=150.):
		def drawObj_multiband(self, tel, pos_img, psf_model, bandpass_list, filt, chip, nphotons_tot=None, g1=0, g2=0,
		exptime=150.):
		if nphotons_tot == None:
		nphotons_tot = self.getElectronFluxFilt(filt, tel, exptime)
		# print("nphotons_tot = ", nphotons_tot)
		@@ -162,7 +164,8 @@ class MockObject(object):
		folding_threshold = 5.e-3
		gsp = galsim.GSParams(folding_threshold=folding_threshold)

		self.real_pos = self.getRealPos(chip.img, global_x=self.posImg.x, global_y=self.posImg.y, img_real_wcs=self.real_wcs)
		self.real_pos = self.getRealPos(chip.img, global_x=self.posImg.x, global_y=self.posImg.y,
		img_real_wcs=self.real_wcs)

		x, y = self.real_pos.x + 0.5, self.real_pos.y + 0.5
		x_nominal = int(np.floor(x + 0.5))
		@@ -192,7 +195,8 @@ class MockObject(object):
		continue
		nphotons_sum += nphotons
		# print("nphotons_sub-band_%d = %.2f"%(i, nphotons))
		psf, pos_shear = psf_model.get_PSF(chip=chip, pos_img=pos_img, bandpass=bandpass, folding_threshold=folding_threshold)
		psf, pos_shear = psf_model.get_PSF(chip=chip, pos_img=pos_img, bandpass=bandpass,
		folding_threshold=folding_threshold)
		star = galsim.DeltaFunction(gsparams=gsp)
		star = star.withFlux(nphotons)
		star = galsim.Convolve(psf, star)
		@@ -205,8 +209,6 @@ class MockObject(object):
		# photons.y += self.y_nominal
		# photons_list.append(photons)



		stamp = star.drawImage(wcs=real_wcs_local, method='phot', offset=offset, save_photons=True)
		xmax = max(xmax, stamp.xmax)
		ymax = max(ymax, stamp.ymax)
		@@ -251,7 +253,6 @@ class MockObject(object):
		del stamp
		return True, pos_shear


		def addSLStoChipImage(self, sdp=None, chip=None, xOrderSigPlus=None, local_wcs=None):
		spec_orders = sdp.compute_spec_orders()
		for k, v in spec_orders.items():
		@@ -322,7 +323,8 @@ class MockObject(object):
		normalSED = Table(np.array([self.sed['WAVELENGTH'], self.sed['FLUX'] * sedNormFactor]).T,
		names=('WAVELENGTH', 'FLUX'))

		self.real_pos = self.getRealPos(chip.img, global_x=self.posImg.x, global_y=self.posImg.y, img_real_wcs=self.real_wcs)
		self.real_pos = self.getRealPos(chip.img, global_x=self.posImg.x, global_y=self.posImg.y,
		img_real_wcs=self.real_wcs)

		x, y = self.real_pos.x + 0.5, self.real_pos.y + 0.5
		x_nominal = int(np.floor(x + 0.5))
		@@ -333,12 +335,15 @@ class MockObject(object):

		real_wcs_local = self.real_wcs.local(self.real_pos)

		flat_cube = chip.flat_cube

		xOrderSigPlus = {'A':1.3909419820029296,'B':1.4760376591236062,'C':4.035447379743442,'D':5.5684364343742825,'E':16.260021029735388}
		xOrderSigPlus = {'A': 1.3909419820029296, 'B': 1.4760376591236062, 'C': 4.035447379743442,
		'D': 5.5684364343742825, 'E': 16.260021029735388}
		grating_split_pos_chip = 0 + grating_split_pos
		for i in range(len(bandpass_list)):
		bandpass = bandpass_list[i]
		psf, pos_shear = psf_model.get_PSF(chip=chip, pos_img=pos_img, bandpass=bandpass, folding_threshold=folding_threshold)
		psf, pos_shear = psf_model.get_PSF(chip=chip, pos_img=pos_img, bandpass=bandpass,
		folding_threshold=folding_threshold)
		star = galsim.DeltaFunction(gsparams=gsp)
		star = star.withFlux(tel.pupil_area * exptime)
		star = galsim.Convolve(psf, star)
		@@ -365,7 +370,8 @@ class MockObject(object):
		tar_spec=normalSED,
		band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10,
		conf=chip.sls_conf[0],
		isAlongY=0)
		isAlongY=0,
		flat_cube=flat_cube)

		self.addSLStoChipImage(sdp=sdp_p1, chip=chip, xOrderSigPlus=xOrderSigPlus, local_wcs=real_wcs_local)

		@@ -380,7 +386,8 @@ class MockObject(object):
		tar_spec=normalSED,
		band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10,
		conf=chip.sls_conf[1],
		isAlongY=0)
		isAlongY=0,
		flat_cube=flat_cube)

		self.addSLStoChipImage(sdp=sdp_p2, chip=chip, xOrderSigPlus=xOrderSigPlus, local_wcs=real_wcs_local)

		@@ -392,7 +399,8 @@ class MockObject(object):
		tar_spec=normalSED,
		band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10,
		conf=chip.sls_conf[1],
		isAlongY=0)
		isAlongY=0,
		flat_cube=flat_cube)
		self.addSLStoChipImage(sdp=sdp, chip=chip, xOrderSigPlus=xOrderSigPlus, local_wcs=real_wcs_local)
		del sdp
		elif grating_split_pos_chip >= gal_end[1]:
		@@ -401,7 +409,8 @@ class MockObject(object):
		tar_spec=normalSED,
		band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10,
		conf=chip.sls_conf[0],
		isAlongY=0)
		isAlongY=0,
		flat_cube=flat_cube)
		self.addSLStoChipImage(sdp=sdp, chip=chip, xOrderSigPlus=xOrderSigPlus, local_wcs=real_wcs_local)
		del sdp
		del psf

ObservationSim/MockObject/SkybackgroundMap.py

+10 −4

Original line number	Diff line number	Diff line
		@@ -20,7 +20,7 @@ except ImportError:
		###calculate sky map by sky SED

		def calculateSkyMap_split_g(skyMap=None, blueLimit=4200, redLimit=6500, skyfn='sky_emiss_hubble_50_50_A.dat', conf=[''], pixelSize=0.074, isAlongY=0,
		split_pos=3685):
		split_pos=3685, flat_cube = None):
		# skyMap = np.ones([yLen, xLen], dtype='float32')
		#
		# if isAlongY == 1:
		@@ -58,6 +58,7 @@ def calculateSkyMap_split_g(skyMap=None, blueLimit=4200, redLimit=6500, skyfn='s
		if split_pos >= skyImg.array.shape[directParm]:
		skyImg1 = galsim.Image(skyImg.array)
		origin1 = [0, 0]
		skyImg1.setOrigin(origin1)
		# sdp = specDisperser.specDisperser(orig_img=skyImg1, xcenter=skyImg1.center.x, ycenter=skyImg1.center.y,
		# full_img=fimg, tar_spec=spec, band_start=tbstart, band_end=tbend,
		# origin=origin1,
		@@ -68,7 +69,8 @@ def calculateSkyMap_split_g(skyMap=None, blueLimit=4200, redLimit=6500, skyfn='s
		sdp = SpecDisperser(orig_img=skyImg1, xcenter=skyImg1.center.x, ycenter=skyImg1.center.y, origin=origin1,
		tar_spec=spec,
		band_start=tbstart, band_end=tbend,
		conf=conf2)
		conf=conf2,
		flat_cube=flat_cube)

		spec_orders = sdp.compute_spec_orders()

		@@ -89,8 +91,10 @@ def calculateSkyMap_split_g(skyMap=None, blueLimit=4200, redLimit=6500, skyfn='s

		skyImg1 = galsim.Image(skyImg.array[:, 0:split_pos])
		origin1 = [0, 0]
		skyImg1.setOrigin(origin1)
		skyImg2 = galsim.Image(skyImg.array[:, split_pos:-1])
		origin2 = [0, split_pos]
		skyImg2.setOrigin(origin2)

		# sdp = specDisperser.specDisperser(orig_img=skyImg1, xcenter=skyImg1.center.x, ycenter=skyImg1.center.y,
		# full_img=fimg, tar_spec=spec, band_start=tbstart, band_end=tbend,
		@@ -102,7 +106,8 @@ def calculateSkyMap_split_g(skyMap=None, blueLimit=4200, redLimit=6500, skyfn='s
		sdp = SpecDisperser(orig_img=skyImg1, xcenter=skyImg1.center.x, ycenter=skyImg1.center.y, origin=origin1,
		tar_spec=spec,
		band_start=tbstart, band_end=tbend,
		conf=conf1)
		conf=conf1,
		flat_cube=flat_cube)

		spec_orders = sdp.compute_spec_orders()

		@@ -121,7 +126,8 @@ def calculateSkyMap_split_g(skyMap=None, blueLimit=4200, redLimit=6500, skyfn='s
		sdp = SpecDisperser(orig_img=skyImg2, xcenter=skyImg2.center.x, ycenter=skyImg2.center.y, origin=origin2,
		tar_spec=spec,
		band_start=tbstart, band_end=tbend,
		conf=conf2)
		conf=conf2,
		flat_cube=flat_cube)

		spec_orders = sdp.compute_spec_orders()