Changeset 13120 for fact/tools/pyscripts/examples
- Timestamp:
- 03/15/12 13:57:33 (13 years ago)
- Location:
- fact/tools/pyscripts/examples
- Files:
-
- 1 added
- 1 deleted
- 2 moved
Legend:
- Unmodified
- Added
- Removed
-
fact/tools/pyscripts/examples/introductory_course/produce_npz.py
r13104 r13120 4 4 # TU Dortmund 5 5 # 6 # test for freshly implemented methods in pyfact.7 # subject to frequent change!6 # example for storing numpy arrays containing FACT analysis results 7 # to file 8 8 9 from pyfact import * 10 #from savecam import * 9 #import this at first, please 10 from pyfact_rename import RawData 11 11 12 import os.path 13 import numpy as np 12 14 13 dfname = '/data00/fact-construction/raw/2011/11/24/20111124_081.fits.gz' 14 outfname = '20111124_081.npz' 15 calfname = '/data00/fact-construction/raw/2011/11/24/20111124_071.drs.fits.gz' 15 from drs_spikes import DRSSpikes 16 from fir_filter import SlidingAverage 17 from extractor import GlobalMaxFinder 18 from cleaners import AmplitudeCleaner 16 19 17 # access the data18 rd = rawdata( dfname, calfname ) 20 ############################################################################## 21 confirm_next_step = False# this is for user interaction 19 22 20 print 'dfname: ', dfname 21 print 'calfname:', calfname 22 print 'NEvents: ', rd.NEvents 23 data_file_name = '/media/daten_platte/FACT/data/20120229_144.fits.gz' 24 if not os.path.isfile(data_file_name): 25 print 'not able to find file:', data_file_name 26 sys.exit(-1) 27 calib_file_name = '/media/daten_platte/FACT/data/20120229_132.drs.fits.gz' 28 if not os.path.isfile(calib_file_name ): 29 print 'not able to find file:', calib_file_name 30 sys.exit(-1) 23 31 24 def loop_acal( Neve = 1000 ): 25 """ 26 bla 27 """ 28 if rd.NEvents < Neve: 29 print 'data file contains not enough events' 30 exit( 0 ) 31 else: 32 for i in range( Neve ): 33 if (np.mod(i,100)==0) : 34 print 'Event: ', i 35 rd.next() 36 37 rd.filterSlidingAverage() 38 # print rd.smoothData 39 rd.filterCFD() 40 # print rd.cfdData 41 rd.findPeak() 42 list_maxPos.append(rd.maxPos) 43 list_maxAmp.append(rd.maxAmp) 44 rd.sumAroundPeak() 45 list_integral.append(rd.integral) 46 # plotincam(rd.integral, 'evt_'+str(i)+'.pdf') 47 #print len(rd.maxPos) 48 #print len(rd.maxAmp) 49 # print 'Trigger Type', rd.trigType 32 run = RawData(data_file_name, calib_file_name) 33 despike = DRSSpikes() 34 smooth = SlidingAverage(8) 35 extract = GlobalMaxFinder(40,200) 36 cleaner = AmplitudeCleaner(45, 18) 50 37 51 #print (rd.filterSlidingAverage.__doc__) 52 #print (rd.filterCFD.__doc__) 53 #print (rd.findPeak.__doc__) 54 #print (rd.sumAroundPeak.__doc__) 38 areas = [] 39 sizes = [] 40 for data,startcell,tt in run: 41 # trigger type 4 means 'physics event' 42 if tt==4: 43 data = despike(data) 44 data = smooth(data) 45 amplitude, time_of_max = extract(data) 46 survivors = cleaner(amplitude, return_bool_mask=False) 55 47 56 list_maxPos = [] 57 list_maxAmp = [] 58 list_integral = [] 59 loop_acal(rd.NEvents) 60 maxAmp = np.vstack(list_maxAmp) 61 maxPos = np.vstack(list_maxPos) 62 integ = np.vstack(list_integral) 48 # this is not python like, should be done in a single line 49 # adding up the amplitude of all survivors 50 size = 0 51 for pixel in survivors: 52 size += amplitude[pixel] 53 54 area = len(survivors) 55 56 if area > 0: 57 areas.append( area ) 58 sizes.append( size ) 59 60 # This is for ---------- USER INTERACTION ------------------------- 61 if confirm_next_step: 62 user_input = raw_input("'q'-quit, 'r'-run, anything else goes one step") 63 number=None 64 try: 65 number=int(user_input) 66 except: 67 number=None 68 if user_input.find('q') != -1: 69 sys.exit(0) 70 elif user_input.find('r') != -1: 71 confirm_next_step = False 72 elif number!=None: 73 run += number 74 # ---------------END OF USER INTERACTION ------------------------- 63 75 64 np.savez(outfname, amplitude=maxAmp, time=maxPos, integral=integ) 76 outfname = raw_input("please enter name for outputfile") 77 np.savez(outfname, Area=areas, Size=sizes)
Note:
See TracChangeset
for help on using the changeset viewer.