Commit 166b2769 authored by Yan Zhaojun's avatar Yan Zhaojun
Browse files

debug

parent 9a787400
Loading
Loading
Loading
Loading
Loading
+7 −167
Original line number Diff line number Diff line
@@ -42,7 +42,6 @@ class cosmicrays():
        :param information: cosmic ray track information (file containing track length and energy information) and
                            exposure time.
        """
        #setup logger
               
        self.exptime=exptime
        
@@ -52,31 +51,13 @@ class cosmicrays():
        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):
@@ -198,113 +179,7 @@ class cosmicrays():
        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):
@@ -361,42 +236,10 @@ class cosmicrays():
            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
            # text = 'The cosmic ray covering factor is %i pixels i.e. %.3f per cent' % (area_cr, covering)
            # self.log.info(text)

###################################################33

    def addUpToFraction(self, coveringFraction, limit=None, verbose=False):
        """
@@ -420,8 +263,5 @@ class cosmicrays():
        return self.image


# if __name__ == "__main__":
    

#     print()