Changeset 13298 for fact/tools
- Timestamp:
- 04/03/12 11:46:51 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
fact/tools/pyscripts/examples/hillas_test.py
r13290 r13298 4 4 # * looping over RawData class object 5 5 # 6 #import matplotlib.pyplot as plt6 import matplotlib.pyplot as plt 7 7 import numpy as np 8 #import matplotlib.patches as patches8 import matplotlib.patches as patches 9 9 10 10 from pyfact import RawData 11 11 #from plotters import Plotter 12 #from plotters import CamPlotter12 from plotters import CamPlotter 13 13 14 14 … … 32 32 clean = AmplitudeCleaner(45,18) 33 33 clean.return_bool_mask = False 34 hillas = HillasParameter( )34 hillas = HillasParameter( source_pos=(1,-1) ) 35 35 36 #p = CamPlotter('cleaned') 37 #p2 = CamPlotter('not') 36 p = CamPlotter('cleaned') 38 37 39 38 39 print 'Source Position for HillasParameter class is:' 40 print hillas.source_pos 41 print 'change it in line 34 if you like' 42 40 43 for event in run: 41 #if event['trigger_type'].value == 4: 42 if True: 43 print event['event_id'] 44 # trigger_type 4 means 'physical data' 45 if event['trigger_type'].value == 4: 46 47 print 'event number:', int(event['event_id'].value) 48 49 #remark: 50 # currently one needs a quite long expression for getting the event id 51 # out of the event-dictionary: 52 # it looks like this: int(event['event_id'].value 53 # 54 # I will change it, such that it looks like this: 55 # event['event_id'] 56 # or maybe: 57 # event['id'] 58 # 59 # what do you think, dear reader?? :-) 60 44 61 data = event['acal_data'] 45 62 unspiked_data = despike(data) … … 48 65 survivors= clean(amp) 49 66 67 # the cleaner also sorts all pixels into groups of islands 68 # clean.islands is a list of lists of pixel_ids 69 # the number of lists is simply the number of islands... 70 num_islands = len(clean.islands) 71 50 72 # if nothing survived the cleaning, just go on 51 if len(clean.islands)== 0:73 if num_islands == 0: 52 74 continue 53 75 54 # play with the if statements here, to look at only those events you would like to analyse 55 #if num_islands >= 2 and len(survivors) > 30: 56 #if num_islands == 1 and len(survivors) > 20: 57 if True: 58 59 print 'calling HillasParameter' 76 77 if num_islands == 1: 78 60 79 hillaspar = hillas(survivors, amp) 61 for k in hillaspar.keys(): 62 if k[0] != '_': 63 print k, hillaspar[k] 80 81 print 'size', hillaspar['size'] 82 print 'alpha', hillaspar['alpha'] 83 print 'delta', hillaspar['delta'] 84 print 'width', hillaspar['width'] 85 print 'length', hillaspar['length'] 86 print 'COG in euclidic coordinates:' 87 print hillaspar['cog_euc'] 88 89 90 p ( amp, survivors ) 91 plt.figure( p.fig_id ) 92 93 # paint arrow to COG 94 plt.gca().add_patch( 95 patches.Polygon( 96 [ 97 (hillaspar['source_pos'][0], hillaspar['source_pos'][1]), 98 (hillaspar['cog_euc'][0], hillaspar['cog_euc'][1]) 99 ] ) ) 100 # paint copy of x-axis through COG 101 plt.axhline( hillaspar['cog_euc'][1] ) 102 103 plt.gca().add_patch( 104 patches.Ellipse( 105 ( hillaspar['cog_euc'][0], hillaspar['cog_euc'][1] ), 106 2*hillaspar['length'], 107 2*hillaspar['width'], 108 hillaspar['delta'] /np.pi * 180, 109 facecolor='none' ) ) 64 110 65 #print 66 #print 'delta:', hillaspar['delta']/np.pi * 180 67 #print 'COG:', hillaspar['cog_euc'] 68 #print 'Mxx/size:', hillaspar['Mxx']/hillaspar['size'] 69 #print 'Mxy/size:', hillaspar['Mxy']/hillaspar['size'] 70 #print 'Myy/size:', hillaspar['Myy']/hillaspar['size'] 111 plt.gca().add_patch( 112 patches.Circle( 113 (hillaspar['source_pos'][0], hillaspar['source_pos'][1] ), 114 0.5 ) ) 71 115 72 73 #p ( amp, survivors ) 74 #plt.figure( p.fig_id ) 75 76 # paint arroy to COG 77 #plt.gca().add_patch( 78 # patches.Polygon( 79 # [ 80 # (hillaspar['source_pos'][0], hillaspar['source_pos'][1]), 81 # (hillaspar['cog_euc'][0], hillaspar['cog_euc'][1]) 82 # ] ) ) 83 # paint copy of x-axis through COG 84 #plt.axhline( hillaspar['cog_euc'][1] ) 85 86 #plt.gca().add_patch( 87 # patches.Ellipse( 88 # ( hillaspar['cog_euc'][0], hillaspar['cog_euc'][1] ), 89 # 2*hillaspar['length'], 90 # 2*hillaspar['width'], 91 # hillaspar['delta'] /np.pi * 180, 92 # facecolor='none' ) ) 93 94 #plt.gca().add_patch( 95 # patches.Circle( 96 # (hillaspar['source_pos'][0], hillaspar['source_pos'][1] ), 97 # 0.5 ) ) 98 99 100 101 102 #p2 (amp) 103 104 #sys.exit(0) 105 106 #answer = raw_input('hit <Enter> to go on .... hit "q" to quit') 107 #if 'q' in answer: 108 # break 116 answer = raw_input('hit <Enter> to go on .... hit "q" to quit') 117 if 'q' in answer: 118 break 109 119 print 'good bye'
Note:
See TracChangeset
for help on using the changeset viewer.