From a5489f19152e36e43a762d183e1fbbe350015473 Mon Sep 17 00:00:00 2001 From: yan Date: Sat, 11 May 2024 11:46:37 +0800 Subject: [PATCH] debug --- csst_ifs_sim/support/IFSinstrumentModel.py | 55 --- csst_ifs_sim/support/__init__.py | 0 csst_ifs_sim/support/cosmicrays.py | 424 --------------------- csst_ifs_sim/support/logger.py | 40 -- 4 files changed, 519 deletions(-) delete mode 100644 csst_ifs_sim/support/IFSinstrumentModel.py delete mode 100644 csst_ifs_sim/support/__init__.py delete mode 100644 csst_ifs_sim/support/cosmicrays.py delete mode 100644 csst_ifs_sim/support/logger.py diff --git a/csst_ifs_sim/support/IFSinstrumentModel.py b/csst_ifs_sim/support/IFSinstrumentModel.py deleted file mode 100644 index c47c72e..0000000 --- a/csst_ifs_sim/support/IFSinstrumentModel.py +++ /dev/null @@ -1,55 +0,0 @@ -""" -VIS Instrument Model -==================== - -The file provides a function that returns VIS related information such as pixel -size, dark current, gain, zeropoint, and sky background. - -:requires: NumPy -:requires: numexpr - -:version: 0.7 -""" - -# import datetime, math -# import numpy as np -# import numexpr as ne - - -def IFSinformation(): - """ - Returns a dictionary describing VIS. The following information is provided (id: value - reference):: - - - - :return: instrument model parameters - :rtype: dict - """ - - ######################################################################################################### - #out = dict(readnoise=4, pixel_size=0.1, dark=0.0008333, fullwellcapacity=90000, bluesize=4000, redsize=6000, readtime=300.) - out=dict() - out.update({'dob' : 0, 'rdose' : 8.0e9, - 'parallelTrapfile' : 'cdm_euclid_parallel.dat', 'serialTrapfile' : 'cdm_euclid_serial.dat', - 'beta_s' : 0.6, 'beta_p': 0.6, 'fwc' : 90000, 'vth' : 1.168e7, 't' : 20.48e-3, 'vg' : 6.e-11, - 'st' : 5.0e-6, 'sfwc' : 730000., 'svg' : 1.0e-10}) - return out - - -def CCDnonLinearityModel(data, beta=6e-7): - """ - - The non-linearity is modelled based on the results presented. - :param data: data to which the non-linearity model is being applied to - :type data: ndarray - - :return: input data after conversion with the non-linearity model - :rtype: float or ndarray - """ - out = data-beta*data**2 - - return out - - -# if __name__ == '__main__': -# print('IFSinstrument') diff --git a/csst_ifs_sim/support/__init__.py b/csst_ifs_sim/support/__init__.py deleted file mode 100644 index e69de29..0000000 diff --git a/csst_ifs_sim/support/cosmicrays.py b/csst_ifs_sim/support/cosmicrays.py deleted file mode 100644 index 632649b..0000000 --- a/csst_ifs_sim/support/cosmicrays.py +++ /dev/null @@ -1,424 +0,0 @@ -""" -Cosmic Rays -=========== - -This simple class can be used to include cosmic ray events to an image. -By default the cosmic ray events are drawn from distributions describing -the length and energy of the events. Such distributions can be generated -for example using Stardust code (http://ieeexplore.ieee.org/stamp/stamp.jsp?arnumber=04636917). -The energy of the cosmic ray events can also be set to constant for -testing purposes. The class can be used to draw a single cosmic ray -event or up to a covering fraction. - -:requires: NumPy -:requires: SciPy - -:version: 0.2 - -""" -import numpy as np -from scipy.interpolate import InterpolatedUnivariateSpline - - -class cosmicrays(): - """ - Cosmic ray generation class. Can either draw events from distributions or - set the energy of the events to a constant. - - :param log: logger instance - :param image: image to which cosmic rays are added to (a copy is made not to change the original numpy array) - :param crInfo: column information (cosmic ray file) - :param information: cosmic ray track information (file containing track length and energy information) and - exposure time. - """ - def __init__(self, log, image, crInfo=None, information=None): - """ - Cosmic ray generation class. Can either draw events from distributions or - set the energy of the events to a constant. - - :param log: logger instance - :param image: image to which cosmic rays are added to (a copy is made not to change the original numpy array) - :param crInfo: column information (cosmic ray file) - :param information: cosmic ray track information (file containing track length and energy information) and - exposure time. - """ - #setup logger - self.log = log - - #image and size - self.image = image.copy() - self.ysize, self.xsize = self.image.shape - - #set up the information dictionary, first with defaults and then overwrite with inputs if given - self.information = (dict(cosmicraylengths='/home/yan/csst-master/data/cdf_cr_length.dat', - cosmicraydistance='/home/yan/csst-master/data/cdf_cr_total.dat', - exptime=565)) - if information is not None: - self.information.update(information) - - if crInfo is not None: - self.cr = crInfo - else: - self._readCosmicrayInformation() - - - # def _readCosmicrayInformation(self): - # self.log.info('Reading in cosmic ray information from %s and %s' % (self.information['cosmicraylengths'], - # self.information['cosmicraydistance'])) - # #read in the information from the files - # crLengths = np.loadtxt(self.information['cosmicraylengths']) - # crDists = np.loadtxt(self.information['cosmicraydistance']) - - # #set up the cosmic ray information dictionary - # self.cr = dict(cr_u=crLengths[:, 0], cr_cdf=crLengths[:, 1], cr_cdfn=np.shape(crLengths)[0], - # cr_v=crDists[:, 0], cr_cde=crDists[:, 1], cr_cden=np.shape(crDists)[0]) - - # return self.cr - - - def _cosmicRayIntercepts(self, lum, x0, y0, l, phi): - """ - Derive cosmic ray streak intercept points. - - :param lum: luminosities of the cosmic ray tracks - :param x0: central positions of the cosmic ray tracks in x-direction - :param y0: central positions of the cosmic ray tracks in y-direction - :param l: lengths of the cosmic ray tracks - :param phi: orientation angles of the cosmic ray tracks - - :return: cosmic ray map (image) - :rtype: nd-array - """ - #create empty array - crImage = np.zeros((self.ysize, self.xsize), dtype=np.float64) - - #x and y shifts - dx = l * np.cos(phi) / 2. - dy = l * np.sin(phi) / 2. - mskdx = np.abs(dx) < 1e-8 - mskdy = np.abs(dy) < 1e-8 - dx[mskdx] = 0. - dy[mskdy] = 0. - - #pixels in x-direction - ilo = np.round(x0.copy() - dx) - msk = ilo < 0. - ilo[msk] = 0 - ilo = ilo.astype(int) - - ihi = 1 + np.round(x0.copy() + dx) - msk = ihi > self.xsize - ihi[msk] = self.xsize - ihi = ihi.astype(int) - - #pixels in y-directions - jlo = np.round(y0.copy() - dy) - msk = jlo < 0. - jlo[msk] = 0 - jlo = jlo.astype(int) - - jhi = 1 + np.round(y0.copy() + dy) - msk = jhi > self.ysize - jhi[msk] = self.ysize - jhi = jhi.astype(int) - - #loop over the individual events - for i, luminosity in enumerate(lum): - n = 0 # count the intercepts - - u = [] - x = [] - y = [] - - #Compute X intercepts on the pixel grid - if ilo[i] < ihi[i]: - for xcoord in range(ilo[i], ihi[i]): - ok = (xcoord - x0[i]) / dx[i] - if np.abs(ok) <= 0.5: - n += 1 - u.append(ok) - x.append(xcoord) - y.append(y0[i] + ok * dy[i]) - else: - for xcoord in range(ihi[i], ilo[i]): - ok = (xcoord - x0[i]) / dx[i] - if np.abs(ok) <= 0.5: - n += 1 - u.append(ok) - x.append(xcoord) - y.append(y0[i] + ok * dy[i]) - - #Compute Y intercepts on the pixel grid - if jlo[i] < jhi[i]: - for ycoord in range(jlo[i], jhi[i]): - ok = (ycoord - y0[i]) / dy[i] - if np.abs(ok) <= 0.5: - n += 1 - u.append(ok) - x.append(x0[i] + ok * dx[i]) - y.append(ycoord) - else: - for ycoord in range(jhi[i], jlo[i]): - ok = (ycoord - y0[i]) / dy[i] - if np.abs(ok) <= 0.5: - n += 1 - u.append(ok) - x.append(x0[i] + ok * dx[i]) - y.append(ycoord) - - #check if no intercepts were found - if n < 1: - xc = int(np.floor(x0[i])) - yc = int(np.floor(y0[i])) - crImage[yc, xc] += luminosity - - #Find the arguments that sort the intersections along the track - u = np.asarray(u) - x = np.asarray(x) - y = np.asarray(y) - - args = np.argsort(u) - - u = u[args] - x = x[args] - y = y[args] - - #Decide which cell each interval traverses, and the path length - for i in range(1, n - 1): - w = u[i + 1] - u[i] - cx = int(1 + np.floor((x[i + 1] + x[i]) / 2.)) - cy = int(1 + np.floor((y[i + 1] + y[i]) / 2.)) - - if 0 <= cx < self.xsize and 0 <= cy < self.ysize: - crImage[cy, cx] += (w * luminosity) - - return crImage - - - # def _drawCosmicRays(self, limit=None): - # """ - # Add cosmic rays to the arrays based on a power-law intensity distribution for tracks. - - # Cosmic ray properties (such as location and angle) are chosen from random Uniform distribution. - # """ - # #estimate the number of cosmics - # cr_n = self.xsize * self.ysize * 0.014 / 43.263316 * 2. - # #scale with exposure time, the above numbers are for the nominal 565s exposure - # cr_n *= (self.information['exptime'] / 565.0) - - # #assume a power-law intensity distribution for tracks - # fit = dict(cr_lo=1.0e3, cr_hi=1.0e5, cr_q=2.0e0) - # fit['q1'] = 1.0e0 - fit['cr_q'] - # fit['en1'] = fit['cr_lo'] ** fit['q1'] - # fit['en2'] = fit['cr_hi'] ** fit['q1'] - - # #pseudo-random numbers taken from a uniform distribution between 0 and 1 - # np.random.seed() - # luck = np.random.rand(int(np.floor(cr_n))) - - # #draw the length of the tracks - # if self.cr['cr_cdfn'] > 1: - # ius = InterpolatedUnivariateSpline(self.cr['cr_cdf'], self.cr['cr_u']) - # self.cr['cr_l'] = ius(luck) - # else: - # self.cr['cr_l'] = np.sqrt(1.0 - luck ** 2) / luck - - # #draw the energy of the tracks - # if self.cr['cr_cden'] > 1: - # ius = InterpolatedUnivariateSpline(self.cr['cr_cde'], self.cr['cr_v']) - # self.cr['cr_e'] = ius(luck) - # else: - # np.random.seed() - # self.cr['cr_e'] = (fit['en1'] + (fit['en2'] - fit['en1']) * - # np.random.rand(int(np.floor(cr_n)))) ** (1.0 / fit['q1']) - - # #Choose the properties such as positions and an angle from a random Uniform dist - # np.random.seed() - # cr_x = self.xsize * np.random.rand(int(np.floor(cr_n))) - - # np.random.seed() - # cr_y = self.ysize * np.random.rand(int(np.floor(cr_n))) - - # np.random.seed() - # cr_phi = np.pi * np.random.rand(int(np.floor(cr_n))) - - # #find the intercepts - # if limit is None: - # self.cosmicrayMap = self._cosmicRayIntercepts(self.cr['cr_e'], cr_x, cr_y, self.cr['cr_l'], cr_phi) - # print ('Number of cosmic ray events:', len(self.cr['cr_e'])) - # else: - # #limit to electron levels < limit - # msk = self.cr['cr_e'] < limit - # print ('Number of cosmic ray events: %i / %i' % (len(self.cr['cr_e'][msk]), int(np.floor(cr_n)))) - # self.cosmicrayMap = self._cosmicRayIntercepts(self.cr['cr_e'][msk], cr_x[msk], cr_y[msk], - # self.cr['cr_l'][msk], cr_phi[msk]) - - # #count the covering factor - # area_cr = np.count_nonzero(self.cosmicrayMap) - # text = 'The cosmic ray covering factor is %i pixels i.e. %.3f per cent' \ - # % (area_cr, 100.*area_cr / (self.xsize*self.ysize)) - # self.log.info(text) - # print (text) - - - # def _drawSingleEvent(self, limit=1000, cr_n=1): - # """ - # Generate a single cosmic ray event and include it to a cosmic ray map (self.cosmicrayMap). - - # :param limit: limiting energy for the cosmic ray event - # :type limit: float - # :param cr_n: number of cosmic ray events to include - # :type cr_n: int - - # :return: None - # """ - # #pseudo-random numbers taken from a uniform distribution between 0 and 1 - # np.random.seed() - # luck = np.random.rand(cr_n) - - # #draw the length of the tracks - # ius = InterpolatedUnivariateSpline(self.cr['cr_cdf'], self.cr['cr_u']) - # self.cr['cr_l'] = ius(luck) - - # #set the energy directly to the limit - # self.cr['cr_e'] = np.asarray([limit, ]*cr_n) - - # #Choose the properties such as positions and an angle from a random Uniform dist - # np.random.seed() - # cr_x = self.xsize * np.random.rand(int(np.floor(cr_n))) - - # np.random.seed() - # cr_y = self.ysize * np.random.rand(int(np.floor(cr_n))) - - # np.random.seed() - # cr_phi = np.pi * np.random.rand(int(np.floor(cr_n))) - - # #find the intercepts - # self.cosmicrayMap = self._cosmicRayIntercepts(self.cr['cr_e'], cr_x, cr_y, self.cr['cr_l'], cr_phi) - - # #count the covering factor - # area_cr = np.count_nonzero(self.cosmicrayMap) - # text = 'The cosmic ray covering factor is %i pixels i.e. %.3f per cent' \ - # % (area_cr, 100.*area_cr / (self.xsize*self.ysize)) - # self.log.info(text) - # print( text) - - - def _drawEventsToCoveringFactor(self, coveringFraction=3.0, limit=1000, verbose=False): - """ - Generate cosmic ray events up to a covering fraction and include it to a cosmic ray map (self.cosmicrayMap). - - :param coveringFraction: covering fraction of cosmic rya events in per cent of total number of pixels - :type coveringFraction: float - :param limit: limiting energy for the cosmic ray event [None = draw from distribution] - :type limit: None or float - :param verbose: print out information to stdout - :type verbose: bool - - - :return: None - """ - self.cosmicrayMap = np.zeros((self.ysize, self.xsize)) - - #how many events to draw at once, too large number leads to exceeding the covering fraction - cr_n = int(295 * self.information['exptime'] / 565. * coveringFraction / 1.4) - - covering = 0.0 - - while covering < coveringFraction: - #pseudo-random numbers taken from a uniform distribution between 0 and 1 - np.random.seed() - luck = np.random.rand(cr_n) - - #draw the length of the tracks - ius = InterpolatedUnivariateSpline(self.cr['cr_cdf'], self.cr['cr_u']) - self.cr['cr_l'] = ius(luck) - - if limit is None: - ius = InterpolatedUnivariateSpline(self.cr['cr_cde'], self.cr['cr_v']) - self.cr['cr_e'] = ius(luck) - else: - #set the energy directly to the limit - self.cr['cr_e'] = np.asarray([limit,]) - - #Choose the properties such as positions and an angle from a random Uniform dist - np.random.seed() - cr_x = self.xsize * np.random.rand(int(np.floor(cr_n))) - - np.random.seed() - cr_y = self.ysize * np.random.rand(int(np.floor(cr_n))) - - np.random.seed() - cr_phi = np.pi * np.random.rand(int(np.floor(cr_n))) - - #find the intercepts - self.cosmicrayMap += self._cosmicRayIntercepts(self.cr['cr_e'], cr_x, cr_y, self.cr['cr_l'], cr_phi) - - #count the covering factor - area_cr = np.count_nonzero(self.cosmicrayMap) - covering = 100.*area_cr / (self.xsize*self.ysize) - - # text = 'The cosmic ray covering factor is %i pixels i.e. %.3f per cent' % (area_cr, covering) - # self.log.info(text) - - # if verbose: - # print( text) - - - # def addCosmicRays(self, limit=None): - # """ - # Include cosmic rays to the image given. - - # :return: image with cosmic rays - # :rtype: ndarray - # """ - # self._drawCosmicRays(limit=limit) - - # #paste cosmic rays - # self.image += self.cosmicrayMap - - # return self.image - - - # def addSingleEvent(self, limit=None): - # """ - # Include a single cosmic ray event to the image given. - - # :return: image with cosmic rays - # :rtype: ndarray - # """ - # self._drawSingleEvent(limit=limit) - - # #paste cosmic rays - # self.image += self.cosmicrayMap - - # return self.image - - - def addUpToFraction(self, coveringFraction, limit=None, verbose=False): - """ - Add cosmic ray events up to the covering Fraction. - - :param coveringFraction: covering fraction of cosmic rya events in per cent of total number of pixels - :type coveringFraction: float - :param limit: limiting energy for the cosmic ray event [None = draw from distribution] - :type limit: None or float - :param verbose: print out information to stdout - :type verbose: bool - - :return: image with cosmic rays - :rtype: ndarray - """ - self._drawEventsToCoveringFactor(coveringFraction, limit=limit, verbose=verbose) - - #paste cosmic rays - self.image += self.cosmicrayMap - - return self.image - - -# if __name__ == "__main__": - - -# print() - diff --git a/csst_ifs_sim/support/logger.py b/csst_ifs_sim/support/logger.py deleted file mode 100644 index b09ef82..0000000 --- a/csst_ifs_sim/support/logger.py +++ /dev/null @@ -1,40 +0,0 @@ -""" -These functions can be used for logging information. - -.. Warning:: logger is not multiprocessing safe. - -:version: 0.3 -""" -import logging -import logging.handlers - - -def setUpLogger(log_filename, loggername='logger'): - """ - Sets up a logger. - - :param: log_filename: name of the file to save the log. - :param: loggername: name of the logger - - :return: logger instance - """ - # create logger - logger = logging.getLogger(loggername) - logger.setLevel(logging.DEBUG) - # Add the log message handler to the logger - handler = logging.handlers.RotatingFileHandler(log_filename) - #maxBytes=20, backupCount=5) - # create formatter - formatter = logging.Formatter('%(asctime)s - %(module)s - %(funcName)s - %(levelname)s - %(message)s') - # add formatter to ch - handler.setFormatter(formatter) - # add handler to logger - - if (logger.hasHandlers()): - logger.handlers.clear() - - logger.addHandler(handler) - - return logger - - -- GitLab