Index: fact/tools/pyscripts/pyfact/pyfact.py
===================================================================
--- fact/tools/pyscripts/pyfact/pyfact.py	(revision 12947)
+++ fact/tools/pyscripts/pyfact/pyfact.py	(revision 13124)
@@ -4,5 +4,5 @@
 # ETH Zurich
 #
-from   ctypes import *
+from ctypes import *
 import numpy as np
 from scipy import signal
@@ -14,7 +14,8 @@
 from ROOT import *
 
-class rawdata( object ):
-    """
-    raw data access and calibration
+
+class RawData( object ):
+    """ raw data access and calibration
+    
     - open raw data file and drs calibration file
     - performs amplitude calibration
@@ -23,22 +24,30 @@
       row = number of pixel
       col = length of region of interest
+      
     """
-    # constructor of the classe
-    def __init__( self, dfname,  calfname, bslfname='' ):
-        """
+
+    def __init__(self, data_file_name, calib_file_name, 
+	    	 user_action_calib=lambda acal_data, data, blm, tom, gm, scells, nroi: None,
+	    	 baseline_file_name=''):
+        """ initialize object
+
         open data file and calibration data file
-        get basic information about the data in dfname
+        get basic information about the data in data_file_name
         allocate buffers for data access
 
-        dfname   : fits or fits.gz file containing the data including the path
-        calfname : fits or fits.gz file containing DRS calibration data
-        bslfname : npy file containing the baseline values
-        """
-        self.dfname   = dfname
-        self.calfname = calfname
-        self.bslfname = bslfname
+        data_file_name   : fits or fits.gz file of the data including the path
+        calib_file_name : fits or fits.gz file containing DRS calibration data
+        baseline_file_name : npy file containing the baseline values
+        
+        """
+
+        self.data_file_name = data_file_name
+        self.calib_file_name = calib_file_name
+        self.baseline_file_name = baseline_file_name
+	    
+	self.user_action_calib = user_action_calib
         
         # baseline correction: True / False
-        if len( bslfname ) == 0:
+        if len(baseline_file_name) == 0:
             self.correct_baseline = False
         else:
@@ -47,208 +56,175 @@
         # access data file
         try:
-            df = fits( self.dfname )
+            data_file = fits(self.data_file_name)
         except IOError:
-            print 'problem accessing data file: ', dfname
-            raise # stop ! no data
-        self.df = df
+            print 'problem accessing data file: ', data_file_name
+            raise  # stop ! no data
+        #: data file (fits object)
+        self.data_file = data_file
         
         # get basic information about the data file
-        self.NROI    = df.GetUInt( 'NROI' ) # region of interest (length of DRS pipeline read out)
-        self.NPIX    = df.GetUInt( 'NPIX' ) # number of pixels (should be 1440)
-        self.NEvents = df.GetNumRows()      # find number of events
+        #: region of interest (number of DRS slices read)
+        self.nroi    = data_file.GetUInt('NROI')
+        #: number of pixels (should be 1440)
+        self.npix    = data_file.GetUInt('NPIX')
+        #: number of events in the data run
+        self.nevents = data_file.GetNumRows()
+        
         # allocate the data memories
-        self.evNum = c_ulong()
-        self.trigType = c_ushort()
-        self.Data  = np.zeros( self.NPIX * self.NROI, np.int16 ) 
-        self.startCells = np.zeros( self.NPIX, np.int16 )
+        self.event_id = c_ulong()
+        self.trigger_type = c_ushort()
+        #: 1D array with raw data
+        self.data  = np.zeros( self.npix * self.nroi, np.int16 )
+        #: slice where drs readout started
+        self.start_cells = np.zeros( self.npix, np.int16 )
+        #: time when the FAD was triggered, in some strange units...
+        self.board_times = np.zeros( 40, np.int32 )
+
         # set the pointers to the data++
-        df.SetPtrAddress( 'EventNum', self.evNum )
-        df.SetPtrAddress( 'TriggerType', self.trigType )
-        df.SetPtrAddress( 'StartCellData', self.startCells ) # DRS readout start cell
-        df.SetPtrAddress( 'Data', self.Data ) # this is what you would expect
-        # df.GetNextRow() # access the first event
-        
-        # access calibration file
+        data_file.SetPtrAddress('EventNum', self.event_id)
+        data_file.SetPtrAddress('TriggerType', self.trigger_type)
+        data_file.SetPtrAddress('StartCellData', self.start_cells) 
+        data_file.SetPtrAddress('Data', self.data) 
+        data_file.SetPtrAddress('BoardTime', self.board_times) 
+                
+        # open the calibration file
         try:
-            calf = fits( self.calfname )
+            calib_file = fits(self.calib_file_name)
         except IOError:
-            print 'problem accessing calibration file: ', calfname
+            print 'problem accessing calibration file: ', calib_file_name
             raise
-        self.calf = calf
-        #
-        BaselineMean      = calf.GetN('BaselineMean')
-        GainMean          = calf.GetN('GainMean')
-        TriggerOffsetMean = calf.GetN('TriggerOffsetMean')
-
-        self.blm = np.zeros( BaselineMean, np.float32 )
-        self.gm  = np.zeros( GainMean, np.float32 )
-        self.tom = np.zeros( TriggerOffsetMean, np.float32 )
-
-        self.Nblm = BaselineMean / self.NPIX
-        self.Ngm  = GainMean / self.NPIX
-        self.Ntom  = TriggerOffsetMean / self.NPIX
-
-        calf.SetPtrAddress( 'BaselineMean', self.blm )
-        calf.SetPtrAddress( 'GainMean', self.gm )
-        calf.SetPtrAddress( 'TriggerOffsetMean', self.tom )
-        calf.GetRow(0)
-
-        self.v_bsl = np.zeros( self.NPIX ) # array with baseline values (all ZERO)
-        self.smoothData = None
-        self.maxPos = None
-        self.maxAmp = None
-
-
-    def next( self ):
-        """
-        load the next event from disk and calibrate it
-        """
-        self.df.GetNextRow()
-        self.calibrate_drsAmplitude()
-
-        
-    def calibrate_drsAmplitude( self ):
-        """
-        perform amplitude calibration for the event 
-        """
-        tomV = 2000./4096.
-        acalData = self.Data * tomV # convert into mV
-
-        # reshape arrays: row = pixel, col = drs_slice
-        acalData = np.reshape( acalData, (self.NPIX, self.NROI) )
-        blm = np.reshape( self.blm, (self.NPIX, 1024) )
-        tom = np.reshape( self.tom, (self.NPIX, 1024) )
-        gm  = np.reshape( self.gm,  (self.NPIX, 1024) )
-        
-        # print 'acal Data ', acalData.shape
-        # print 'blm shape ', blm.shape
-        # print 'gm shape  ', gm.shape
-        
-        for pixel in range( self.NPIX ):
+        #: drs calibration file
+        self.calib_file = calib_file
+        
+        baseline_mean = calib_file.GetN('BaselineMean')
+        gain_mean = calib_file.GetN('GainMean')
+        trigger_offset_mean = calib_file.GetN('TriggerOffsetMean')
+
+        self.blm = np.zeros(baseline_mean, np.float32)
+        self.gm  = np.zeros(gain_mean, np.float32)
+        self.tom = np.zeros(trigger_offset_mean, np.float32)
+
+        self.Nblm = baseline_mean / self.npix
+        self.Ngm  = gain_mean / self.npix
+        self.Ntom  = trigger_offset_mean / self.npix
+
+        calib_file.SetPtrAddress('BaselineMean', self.blm)
+        calib_file.SetPtrAddress('GainMean', self.gm)
+        calib_file.SetPtrAddress('TriggerOffsetMean', self.tom)
+        calib_file.GetRow(0)
+
+        self.v_bsl = np.zeros(self.npix)  # array of baseline values (all ZERO)
+        self.data_saverage_out = None
+        self.pulse_time_of_maximum = None
+        self.pulse_amplitude = None
+
+    def __iter__(self):
+        """ iterator """
+        return self
+        
+    def __add__(self, jump_over):
+        self.data_file.GetRow(jump_over)
+        return self
+        
+    def next(self):
+        """ used by __iter__ """
+
+        if self.data_file.GetNextRow() == False:
+            raise StopIteration
+        else:
+            self.calibrate_drs_amplitude()
+
+        #print 'nevents = ', self.nevents, 'event_id = ', self.event_id.value
+	
+        return self.acal_data, self.start_cells, self.trigger_type.value
+
+    def next_event(self):
+        """ load the next event from disk and calibrate it
+        
+        """
+
+        self.data_file.GetNextRow()
+        self.calibrate_drs_amplitude()
+
+    def calibrate_drs_amplitude(self):
+        """ perform the drs amplitude calibration of the event data
+        
+        """
+
+        to_mV = 2000./4096.
+        #: 2D array with amplitude calibrated dat in mV
+        acal_data = self.data * to_mV  # convert ADC counts to mV
+
+        # make 2D arrays: row = pixel, col = drs_slice
+        acal_data = np.reshape(acal_data, (self.npix, self.nroi) )
+        blm = np.reshape(self.blm, (self.npix, self.Nblm) )
+        tom = np.reshape(self.tom, (self.npix, self.Ntom) )
+        gm  = np.reshape(self.gm,  (self.npix, self.Ngm) )
+        
+        for pixel in range( self.npix ):
             # rotate the pixel baseline mean to the Data startCell
-            blm_pixel = np.roll( blm[pixel,:], -self.startCells[pixel] )
-            acalData[pixel,:] -= blm_pixel[0:self.NROI]
-            acalData[pixel,:] -= tom[pixel, 0:self.NROI]
-            acalData[pixel,:] /= gm[pixel,  0:self.NROI]
+            blm_pixel = np.roll( blm[pixel,:], -self.start_cells[pixel] )
+            tom_pixel = np.roll( tom[pixel,:], -self.start_cells[pixel] )
+            gm_pixel = np.roll( gm[pixel,:], -self.start_cells[pixel] )
+            acal_data[pixel,:] -= blm_pixel[0:self.nroi]
+            acal_data[pixel,:] -= tom_pixel[0:self.nroi]
+            acal_data[pixel,:] /= gm_pixel[0:self.nroi]
             
-        self.acalData = acalData * 1907.35
-    
-        # print 'acalData ', self.acalData[0:2,0:20]
-        
-    def filterSlidingAverage( self , windowSize = 4):
-        """ sliding average filter
-        using:
-            self.acalData
-        filling array:
-            self.smoothData
-        """
-        #scipy.signal.lfilter(b, a, x, axis=-1, zi=None)
-        smoothData = self.acalData.copy()
-        b = np.ones( windowSize )
-        a = np.zeros( windowSize )
-        a[0] = len(b)
-        smoothData[:,:] = signal.lfilter(b, a, smoothData[:,:])
-
-        self.smoothData = smoothData
-        
-    def filterCFD( self, length=10, ratio=0.75):
-        """ constant fraction filter
-        using:
-            self.smoothData
-        filling array:
-            self.cfdData
-        """
-        if self.smoothData == None:
-            print 'error pyfact.filterCFD was called without prior call to filterSlidingAverage'
-            print ' variable self.smoothData is needed '
-            pass
-        cfdData = self.smoothData.copy()
-        b = np.zeros( length )
-        a = np.zeros( length )
-        b[0] = -1. * ratio
-        b[length-1] = 1.
-        a[0] = 1.
-        cfdData[:,:] = signal.lfilter(b, a, cfdData[:,:])
-        
-        self.cfdData = cfdData
-    
-    def findPeak (self, min=30, max=250):
-        """ find maximum in search window
-        using: 
-            self.smoothData
-        filling arrays:
-            self.maxPos
-            self.maxAmp
-        """
-        if self.smoothData == None:
-            print 'error pyfact.findPeakMax was called without prior call to filterSlidingAverage'
-            print ' variable self.smoothData is needed '
-            pass
-        maxPos = np.argmax( self.smoothData[:,min:max] , 1)
-        maxAmp = np.max( self.smoothData[:,min:max] , 1)
-        self.maxPos = maxPos
-        self.maxAmp = maxAmp
-
-    def sumAroundPeak (self, left=13, right=23):
-        """ integrate signal in gate around Peak
-        using:
-            self.maxPos
-            self.acalData
-        filling array:
-            self.integral
-        """
-        if self.maxPos == None:
-            print 'error pyfact.sumAroundPeak was called without prior call of findPeak'
-            print ' variable self.maxPos is needed'
-            pass
-        
-        sums = np.empty( self.NPIX )
-        for pixel in range( self.NPIX ):
-            min = self.maxPos[pixel]-left
-            max = self.maxPos[pixel]+right
-            sums[pixel] = self.acalData[pixel,min:max].sum()
-        
-        self.integral = sums
-        
-    def ReadBaseline( self, file, bsl_hist = 'bsl_sum/hplt_mean' ):
-        """
+        self.acal_data = acal_data * 1907.35
+	
+	#print 'blm _pyfact', blm[0,0:20] 
+	#t = np.roll( blm[0,:], -self.start_cells[0] )
+	#print 'blm _pyfact', t[0:20]
+	#print 'start_pyfact: ', self.start_cells[0]
+	#print 'acal _pyfact: ', self.acal_data[0,0:10] 
+	#t = np.roll( gm[0,:], -self.start_cells[0] )
+	#print 'gm _pyfact: ', t[0:10] 
+	self.user_action_calib( self.acal_data, 
+		np.reshape(self.data, (self.npix, self.nroi) ), blm, tom, gm, self.start_cells, self.nroi)
+
+        
+    def baseline_read_values(self, file, bsl_hist='bsl_sum/hplt_mean'):
+        """
+        
         open ROOT file with baseline histogram and read baseline values
         file       name of the root file
         bsl_hist   path to the histogram containing the basline values
-        """
+
+        """
+
         try:
-            f = TFile( file )
+            f = TFile(file)
         except:
             print 'Baseline data file could not be read: ', file
             return
         
-        h = f.Get( bsl_hist )
-
-        for i in range( self.NPIX ):
-            self.v_bsl[i] = h.GetBinContent( i+1 )
+        h = f.Get(bsl_hist)
+
+        for i in range(self.npix):
+            self.v_bsl[i] = h.GetBinContent(i+1)
 
         f.Close()
-
-        
-    def CorrectBaseline( self ):
-        """
-        apply baseline correction
-        """
-        for pixel in range( self.NPIX ):
-            self.acalData[pixel,:] -= self.v_bsl[pixel]
-            
-        
-    def info( self ):
-        """
-        print information
-        """
-        print 'data file:  ', dfname
-        print 'calib file: ', calfname
+        
+    def baseline_correct(self):
+        """ subtract baseline from the data
+
+        """
+        
+        for pixel in range(self.npix):
+            self.acal_data[pixel,:] -= self.v_bsl[pixel]
+                    
+    def info(self):
+        """ print run information
+        
+        """
+        
+        print 'data file:  ', data_file_name
+        print 'calib file: ', calib_file_name
         print 'calibration file'
-        print 'N BaselineMean: ', self.Nblm
-        print 'N GainMean: ', self.Ngm
+        print 'N baseline_mean: ', self.Nblm
+        print 'N gain mean: ', self.Ngm
         print 'N TriggeroffsetMean: ', self.Ntom
         
-# --------------------------------------------------------------------------------
+# -----------------------------------------------------------------------------
 class fnames( object ):
     """ organize file names of a FACT data run
@@ -256,5 +232,5 @@
     """
     
-    def __init__( self, specifier = ['012', '023', '2011', '11', '24'],
+    def __init__(self, specifier = ['012', '023', '2011', '11', '24'],
                  rpath = '/scratch_nfs/res/bsl/',
                  zipped = True):
@@ -265,5 +241,7 @@
         rpath     : directory path for the results; YYYYMMDD will be appended to rpath
         zipped    : use zipped (True) or unzipped (Data) 
-        """
+
+        """
+        
         self.specifier = specifier
         self.rpath     = rpath
@@ -271,5 +249,5 @@
         
         self.make( self.specifier, self.rpath, self.zipped )
-    # end of def __init__
+
 
     def make( self, specifier, rpath, zipped ):
@@ -312,6 +290,4 @@
         self.results = self.names['results']
 
-    # end of make
-
     def info( self ):
         """ print complete filenames
@@ -323,61 +299,22 @@
         print 'drs-cal: ', self.names['drscal']
         print 'results: ', self.names['results']
-    # end of def info
 
 # end of class definition: fnames( object )
 
-
-
-class histogramList( object ):
-
-    def __init__( self, name ):
-        """ set the name and create empty lists """
-        self.name  = name         # name of the list
-        self.list  = []           # list of the histograms
-        self.dict  = {}           # dictionary of histograms
-        self.hList = TObjArray()  # list a la ROOT of the histograms
-
-    def add( self, tag, h ):
-        self.list.append( h )
-        self.dict[tag] = h
-        self.hList.Add( h )
-
-
-class pixelHisto1d ( object ):
-
-    def __init__( self, name, title, Nbin, first, last, xtitle, ytitle, NPIX ):
-        """
-        book one dimensional histograms for each pixel
-        """
-        self.name = name
-
-        self.list = [ x for x in range( NPIX ) ]
-        self.hList = TObjArray()
-
-        for pixel in range( NPIX ):
-
-            hname  = name + ' ' + str( pixel )
-            htitle = title + ' ' + str( pixel )
-            self.list[pixel] = TH1F( hname, htitle, Nbin, first, last )
-
-            self.list[pixel].GetXaxis().SetTitle( xtitle )
-            self.list[pixel].GetYaxis().SetTitle( ytitle )
-            self.hList.Add( self.list[pixel] )
-
-# simple test method
+def _test_iter():
+    """ test for function __iter__ """
+
+#    data_file_name = '/data00/fact-construction/raw/2011/11/24/20111124_111.fits.gz'
+#    calib_file_name =     data_file_name = 
+    data_file_name =  '/home/luster/win7/FACT/data/raw/20120114/20120114_028.fits.gz'
+    calib_file_name = '/home/luster/win7/FACT/data/raw/20120114/20120114_022.drs.fits.gz'
+    run = RawData( data_file_name, calib_file_name )
+
+    for data, scell, tt in run:
+        print 'data[0,0] = ', data[0,0], "start_cell[0] =", scell[0], "trigger type = ", tt
+
+
 if __name__ == '__main__':
-    """
-    create an instance
-    """
-    dfname = '/data03/fact-construction/raw/2011/11/24/20111124_121.fits'
-    calfname = '/data03/fact-construction/raw/2011/11/24/20111124_111.drs.fits'
-    rd = rawdata( dfname, calfname )
-    rd.info()
-    rd.next()
-    
-# for i in range(10):
-#    df.GetNextRow() 
-
-#    print 'evNum: ', evNum.value
-#    print 'startCells[0:9]: ', startCells[0:9]
-#    print 'evData[0:9]: ', evData[0:9]
+    """ tests  """
+
+    _test_iter()
