Chip_deprecated.py

from .FocalPlane import FocalPlane
import galsim
import os
import numpy as np
from astropy.table import Table

class Chip(FocalPlane):
	def __init__(self, chipID, ccdEffCurve_dir, config=None, treering_func=None):
		# Get focal plane (instance of paraent class) info
		# TODO: use chipID to config individual chip?
		super().__init__()
		if config is not None:
			self.npix_x 	= config["npix_x"]
			self.npix_y 	= config["npix_y"]
			self.read_noise = config["read_noise"]
			self.dark_noise = config["dark_noise"]
			self.pix_scale 	= config["pix_scale"]
			self.gain		= config["gain"]
		else:
			# Default setting
			self.npix_x 	= 9232
			self.npix_y 	= 9216
			self.read_noise = 5.0	# e/pix
			self.dark_noise = 0.02	# e/pix/s
			self.pix_scale 	= 0.074	# pixel scale
			self.gain		= 1.0

		# A chip ID must be assigned
		self.chipID = int(chipID)

		# Get corresponding filter info
		self.filter_id, self.filter_type = self.getChipFilter()

		# Get boundary (in pix)
		self.bound = self.getChipLim()
		self.ccdEffCurve_dir = ccdEffCurve_dir
		self.effCurve = self._getChipEffCurve(self.filter_type)

		# Define the sensor
		# self.sensor = galsim.Sensor()
		self.sensor = galsim.SiliconSensor(strength=10., treering_func=treering_func)

	def _getChipEffCurve(self, filter_type):
		if filter_type in ['nuv', 'u']: filename = 'UV0.txt'
		if filter_type in ['g', 'r']: filename = 'Astro_MB.txt'
		if filter_type in ['i', 'z', 'y']: filename = 'Basic_NIR.txt'
		path = os.path.join(self.ccdEffCurve_dir, filename)
		table = Table.read(path, format='ascii')
		throughput = galsim.LookupTable(x=table['col1'], f=table['col2'], interpolant='linear')
		bandpass = galsim.Bandpass(throughput, wave_type='nm')
		return bandpass

	def getChipFilter(self, chipID=None, filter_layout=None):
		"""Return the filter index and type for a given chip #(chipID)
		"""
		filter_type_list = ["nuv","u", "g", "r", "i","z","y"]
		# TODO: maybe a more elegent way other than hard coded?
		# e.g. use something like a nested dict:
		if filter_layout is not None:
			return filter_layout[chipID][0], filter_layout[chipID][1]
		if chipID == None:
			chipID = self.chipID

		# updated configurations
		if chipID>30 or chipID<1: raise ValueError("!!! Chip ID: [1,30]")
		if chipID in [10, 11, 20, 21]: filter_type = 'y'
		if chipID in [15, 16]:         filter_type = "z"
		if chipID in [9, 22]:			filter_type = "i"
		if chipID in [12, 19]:         filter_type = "u"
		if chipID in [7, 24]:         filter_type = "r"
		if chipID in [14, 13, 18, 17]:    filter_type = "nuv"
		if chipID in [8, 23]:         filter_type = "g"
		if chipID in [6, 5, 25, 26]:	filter_type = "GI"
		if chipID in [27, 30, 1, 4]:	filter_type = "GV"
		if chipID in [28, 29, 2, 3]:	filter_type = "GU"
		filter_id = filter_type_list.index(filter_type)

		return filter_id, filter_type

	def getChipLim(self, chipID=None):
		"""Calculate the edges in pixel for a given CCD chip on the focal plane
		NOTE: There are 5*4 CCD chips in the focus plane for photometric observation.
		Parameters:
			chipID:			int
							the index of the chip
		Returns:
			A galsim BoundsD object
		"""
		if chipID == None:
			chipID = self.chipID
		
		gx = self.npix_gap_x
		gy1, gy2 = self.npix_gap_y

		# xlim of a given ccd chip
		xrem = (chipID-1)%self.nchip_x - self.nchip_x // 2
		xcen = (self.npix_x + gx) * xrem
		nx0 = xcen - self.npix_x//2 + 1
		nx1 = xcen + self.npix_x//2

		# ylim of a given ccd chip
		yrem = 2*((chipID-1)//self.nchip_x) - (self.nchip_y-1)
		ycen = (self.npix_y//2 + gy1//2) * yrem
		if chipID <= 6: ycen = (self.npix_y//2 + gy1//2) * yrem - (gy2-gy1)
		if chipID >= 25: ycen = (self.npix_y//2 + gy1//2) * yrem + (gy2-gy1)
		ny0 = ycen - self.npix_y//2 + 1
		ny1 = ycen + self.npix_y//2

		return galsim.BoundsD(nx0-1, nx1-1, ny0-1, ny1-1)


	def getSkyCoverage(self, wcs):
		return super().getSkyCoverage(wcs, self.bound.xmin, self.bound.xmax, self.bound.ymin, self.bound.ymax)


	def getSkyCoverageEnlarged(self, wcs, margin=0.5):
		"""The enlarged sky coverage of the chip
		"""
		margin /= 60.0
		bound = self.getSkyCoverage(wcs)
		return galsim.BoundsD(bound.xmin - margin, bound.xmax + margin, bound.ymin - margin, bound.ymax + margin)


	def isContainObj(self, ra_obj, dec_obj, sky_coverage):
		if (ra_obj - sky_coverage.xmin) * (ra_obj - sky_coverage.xmax) > 0.0 or (dec_obj - sky_coverage.ymin) * (dec_obj - sky_coverage.ymax) > 0.0:
			return False
		return True


	def getChipNoise(self, exptime=150.0):
		# TODO: put it here or make it a separate class?
		noise = self.dark_noise * exptime + self.read_noise**2
		return noise


	def addNoise(self, img, exptime=150.0, sky_noise=0., seed=31415):
		rng = galsim.BaseDeviate(seed)
		dark_noise = galsim.DeviateNoise(galsim.PoissonDeviate(rng, self.dark_noise*exptime))
		img.addNoise(dark_noise)
		ccd_noise = galsim.CCDNoise(rng, sky_level=sky_noise, gain=self.gain, read_noise=self.read_noise)
		img.addNoise(ccd_noise)
		return img