source: fact/tools/pyscripts/new_pyfact/pyfact.py@ 20007

Last change on this file since 20007 was 17979, checked in by dneise, 11 years ago
removed all pylint warnings and errors
  • Property svn:executable set to *
File size: 13.8 KB
Line 
1#!/usr/bin/python -tt
2#
3# Werner Lustermann, Dominik Neise
4# ETH Zurich, TU Dortmund
5#
6import os
7import sys
8import numpy as np
9import ROOT
10ROOT.gROOT.SetBatch(True)
11
12#--------- BUILDING OF THE SHARED OBJECT FILES -------------------------------
13if __name__ == '__main__' and len(sys.argv) > 1 and 'build' in sys.argv[1]:
14 ROOT.gSystem.AddLinkedLibs("-lz")
15 root_make_string = ROOT.gSystem.GetMakeSharedLib()
16 if "-std=c++0x" not in root_make_string:
17 make_string = root_make_string.replace('$Opt', '$Opt -std=c++0x -D HAVE_ZLIB')
18 ROOT.gSystem.SetMakeSharedLib(make_string)
19 ROOT.gROOT.ProcessLine(".L extern_Mars_mcore/izstream.h+O")
20 ROOT.gROOT.ProcessLine(".L extern_Mars_mcore/fits.h+O")
21 ROOT.gROOT.ProcessLine(".L extern_Mars_mcore/zfits.h+O")
22 ROOT.gROOT.ProcessLine(".L extern_Mars_mcore/factfits.h+O")
23 if "-std=c++0x" not in root_make_string:
24 make_string = root_make_string.replace(
25 '$Opt',
26 "$Opt -std=c++0x "
27 "-D HAVE_ZLIB "
28 "-D'PACKAGE_NAME=\"PACKAGE_NAME\"' "
29 "-D'PACKAGE_VERSION=\"PACKAGE_VERSION\"' "
30 "-D'REVISION=\"REVISION\"' "
31 )
32 ROOT.gSystem.SetMakeSharedLib(make_string)
33 ROOT.gROOT.ProcessLine(".L extern_Mars_mcore/DrsCalib.h+O")
34
35 ROOT.gInterpreter.GenerateDictionary(
36 "map<string,fits::Entry>",
37 "map;string;extern_Mars_mcore/fits.h")
38 ROOT.gInterpreter.GenerateDictionary(
39 "pair<string,fits::Entry>",
40 "map;string;extern_Mars_mcore/fits.h")
41 ROOT.gInterpreter.GenerateDictionary(
42 "map<string,fits::Table::Column>",
43 "map;string;extern_Mars_mcore/fits.h")
44 ROOT.gInterpreter.GenerateDictionary(
45 "pair<string,fits::Table::Column>",
46 "map;string;extern_Mars_mcore/fits.h")
47 ROOT.gInterpreter.GenerateDictionary(
48 "vector<DrsCalibrate::Step>",
49 "vector;extern_Mars_mcore/DrsCalib.h")
50
51#-------- USAGE -------------------------------------------------------------
52if __name__ == '__main__' and len(sys.argv) < 3:
53 print """ Usage:
54 ----------------------------------------------------------------------
55 To just build the shared object libs call:
56 python pyfact.py build
57
58 To build (if necessary) and open an example file
59 python -i pyfact.py /path/to/data_file.zfits /path/to/calib_file.drs.fits.gz
60
61 Any of the 3 file 'types': fits.gz zfits fits
62 should be supported.
63
64 To clean all of automatically built files, do something like:
65 rm *.so *.d AutoDict_* extern_Mars_mcore/*.so extern_Mars_mcore/*.d
66 ----------------------------------------------------------------------
67 """
68 sys.exit(1)
69
70#-------- START OF PYFACT MODULE --------------------------------------------
71path = os.path.dirname(os.path.realpath(__file__))
72ROOT.gSystem.Load(path+'/AutoDict_map_string_fits__Entry__cxx.so')
73ROOT.gSystem.Load(path+'/AutoDict_map_string_fits__Table__Column__cxx.so')
74ROOT.gSystem.Load(path+'/AutoDict_pair_string_fits__Entry__cxx.so')
75ROOT.gSystem.Load(path+'/AutoDict_pair_string_fits__Table__Column__cxx.so')
76ROOT.gSystem.Load(path+'/AutoDict_vector_DrsCalibrate__Step__cxx.so')
77ROOT.gSystem.Load(path+'/extern_Mars_mcore/fits_h.so')
78ROOT.gSystem.Load(path+'/extern_Mars_mcore/izstream_h.so')
79ROOT.gSystem.Load(path+'/extern_Mars_mcore/zfits_h.so')
80ROOT.gSystem.Load(path+'/extern_Mars_mcore/factfits_h.so')
81ROOT.gSystem.Load(path+'/extern_Mars_mcore/DrsCalib_h.so')
82del path
83
84env_name = 'FACT_DATA_BASE_DIR'
85if env_name in os.environ:
86 base_path = os.environ[env_name]
87 raw_base = os.path.join(base_path, 'raw')
88 aux_base = os.path.join(base_path, 'aux')
89else:
90 raw_base = os.path.join('fact', 'raw')
91 aux_base = os.path.join('fact', 'aux')
92
93
94def stub_to_auxfile(s, s2=None, basep=aux_base, basen="*"):
95 """ FIXME :-)
96 It's not yet clear how to generate aux file names from stubs, or iterables
97 """
98 if s2 is not None:
99 s = "{0}_{1:03d}".format(s, s2)
100 y = s[0:4]
101 m = s[4:6]
102 d = s[6:8]
103 date = s[0:8]
104 #r = s[9:12]
105 return os.path.join(basep, y, m, d, date + '.' + basen + '.fits')
106
107
108def make_path(s):
109 """ make path to fits file
110
111 s : string or iterable
112
113 Fits files in pyfact do not need to be specified by full paths
114 A so called stub path is sufficient. A valid stub path looks like:
115 yyyymmdd_rrr
116 or in case of a drs calibration file like
117 yyyymmdd_rrr.drs
118
119 In addition a file can be specified by a 2-element integer iterable,
120 e.g. a list or a tuple, whose first element is the integer denoting the
121 'night' : yyyymmdd
122 the 2nd integer is the 'run' : rrr
123 so valid specifiers are also:
124 (yyyymmdd, rrr) or [yyyymmdd, rrr]
125
126 the pyfact.Fits class can handle .fits, .fits.gz, and .fits.fz files.
127 stubs dont's specify which one should be used, so we need to check
128 which one exists, if more than one file exists fitting to the stub,
129 the order of precendence is:
130 * .fits -- since access should be pretty fast
131 * .fits.fz -- newer format should be mostly the case
132 * .fits.gz -- still possible
133 """
134 import re
135 import os
136 if not isinstance(s, (str, unicode)):
137 # s is (hopefully) some kind of iterable, so we read out two integers
138 night = int(s[0])
139 run = int(s[1])
140 _s = "{0:08d}_{1:03d}".format(night, run)
141 else:
142 _s = s
143
144 # so now s is some kind of string
145 # maybe it is already a path:
146 if os.path.isabs(_s):
147 # s is already an absolute path.
148 # maybe the file is not there, but we don't care about that
149 return _s
150
151 # it was not an absolute path, maybe it's a relative path, let's check if
152 # it points to a file.
153 elif os.path.exists(_s):
154 # s is certainly a path, it even points to a file :-)
155 return os.path.abspath(_s)
156
157 # okay ... s was not a path, so let's generate a nice absolute path
158 # from the stub, in case it is a stub ..
159 elif re.match('\d{8}_\d{3}', _s):
160 # s is a stub, so we can make a full path from it
161 y = _s[0:4]
162 m = _s[4:6]
163 d = _s[6:8]
164 #r = _s[9:12]
165 start = os.path.join(raw_base, y, m, d, _s[:12] + '.fits')
166 candidates = [start, start+'.fz', start+'.gz']
167 for c in candidates:
168 if os.path.exists(c):
169 return c
170 # if we reached this point, s was a stub, but we were not
171 # able to find a file on the file system...
172 # let's return our best guess, the '.fz' file path
173 # the Fits.__init__() will complain, if it can't find the file.
174 return candidates[1]
175 else:
176 raise IOError("{0} can't be used for path generation".format(s))
177
178
179class Fits(object):
180 """ General FITS file access
181
182 Wrapper for factfits class from Mars/mcore/factfits.h
183 (factfits might not be suited to read any FITS file out there, but certainly
184 all FACT FITS files can be read using this class)
185 """
186 __module__ = 'pyfact'
187
188 def __init__(self, path):
189 """ Open fits file, parse header and setup table colums
190
191 path : string : path to the fits file (.fits, .fits.gz, .zfits)
192 """
193 path = make_path(path)
194 self._path = path
195 self._current_row = None
196 if not os.path.exists(path):
197 raise IOError(path+' was not found')
198 self.f = ROOT.factfits(path)
199 self._make_header()
200 self._setup_columns()
201
202 def __repr__(self):
203 return '{s.__class__.__name__}(path="{s._path}")'.format(s=self)
204
205 def __str__(self):
206 s = """Fits object:
207 path = {s._path}
208 row = {s._current_row}
209 self.cols.keys = {col_names}
210 """.format(s=self, col_names=self.cols.keys())
211 return s
212
213 def __len__(self):
214 """ return the number of events, in case it's a FACT physics data file
215 or simply the number of rows in case it's a slow data (so called aux) file.
216 """
217 return self.f.GetNumRows()
218
219 def _make_header(self):
220 """
221 """
222 str_to_bool = {'T': True, 'F': False}
223 type_conversion = {'I': int, 'F': float, 'T': str, 'B': str_to_bool.__getitem__}
224
225 self.header = {}
226 self.header_comments = {}
227 for key, entry in self.f.GetKeys():
228 try:
229 self.header[key] = type_conversion[entry.type](entry.value)
230 except KeyError:
231 raise IOError(
232 "Error: entry type unknown.\n "
233 "Is %s, but should be one of: [I,F,T,B]" % (entry.type))
234 self.header_comments[key] = entry.comment
235
236 def _setup_columns(self):
237 """
238 """
239 col_type_to_np_type_map = {
240 'L': 'b1',
241 'A': 'a1',
242 'B': 'i1',
243 'I': 'i2',
244 'J': 'i4',
245 'K': 'i8',
246 'E': 'f4',
247 'D': 'f8'}
248 self.cols = {}
249 for key, col in self.f.GetColumns():
250 self.cols[key] = np.zeros(col.num, col_type_to_np_type_map[col.type])
251 if col.num != 0:
252 self.f.SetPtrAddress(key, self.cols[key])
253
254 def __iter__(self):
255 return self
256
257 def next(self, row=None):
258 """
259 """
260 if row is None:
261 if self._current_row is None:
262 self._current_row = 0
263 else:
264 self._current_row += 1
265 if self.f.GetNextRow() is False:
266 self._current_row = None
267 raise StopIteration
268 else:
269 row = int(row)
270 self._current_row = row
271 if self.f.GetRow(row) is False:
272 self._current_row = None
273 raise StopIteration
274 return self
275
276
277class AuxFile(Fits):
278 """ easy(?) access to FACT aux files
279 """
280 __module__ = 'pyfact'
281
282 def __init__(self, path, verbose=False):
283 self._verbose = verbose
284 super(AuxFile, self).__init__(path)
285
286 def _setup_columns(self):
287 col_type_to_np_type_map = {
288 'L': 'b1',
289 'A': 'a1',
290 'B': 'i1',
291 'I': 'i2',
292 'J': 'i4',
293 'K': 'i8',
294 'E': 'f4',
295 'D': 'f8'}
296 self._cols = {}
297 self.cols = {}
298
299 N = len(self)
300 for key, col in self.f.GetColumns():
301 self._cols[key] = np.zeros(col.num, col_type_to_np_type_map[col.type])
302 self.cols[key] = np.zeros((N, col.num), col_type_to_np_type_map[col.type])
303 if col.num != 0:
304 self.f.SetPtrAddress(key, self._cols[key])
305
306 for i, row in enumerate(self):
307 if self._verbose:
308 try:
309 step = int(self._verbose)
310 except:
311 step = 10
312 if i % step == 0:
313 print "reading line", i, 'of', self.header['NAXIS2']
314
315 for key in self._cols:
316 self.cols[key][i, :] = self._cols[key]
317
318
319class RawData(Fits):
320 """ Special raw data FITS file access (with DRS4 calibration)
321
322 During iteration the C++ method DrsCalibration::Apply is being called.
323 """
324 __module__ = 'pyfact'
325
326 def __init__(self, data_path, calib_path):
327 """ -constructor-
328 *data_path* : fits or fits.gz file of the data including the path
329 *calib_path* : fits or fits.gz file containing DRS calibration data
330 """
331 super(RawData, self).__init__(data_path)
332 self.cols['CalibData'] = np.zeros(self.cols['Data'].shape, np.float32)
333 self.cols['CalibData2D'] = self.cols['CalibData'].reshape(self.header['NPIX'], -1)
334 if not self.cols['CalibData2D'].base is self.cols['CalibData']:
335 print "Error seomthing went wrong!"
336 self.drs_calibration = ROOT.DrsCalibration()
337 self.drs_calibration.ReadFitsImp(calib_path)
338
339 self.drs_calibrate = ROOT.DrsCalibrate()
340 self.list_of_previous_start_cells = []
341
342 def __iter__(self):
343 """ iterator """
344 return self
345
346 def next(self, row=None):
347 """
348 """
349 super(RawData, self).next(row)
350
351 self.drs_calibration.Apply(
352 self.cols['CalibData'],
353 self.cols['Data'],
354 self.cols['StartCellData'],
355 self.header['NROI']
356 )
357
358 for previous_start_cells in self.list_of_previous_start_cells:
359 self.drs_calibrate.CorrectStep(
360 self.cols['CalibData'],
361 self.header['NPIX'],
362 self.header['NROI'],
363 previous_start_cells,
364 self.cols['StartCellData'],
365 self.header['NROI']+10)
366 self.drs_calibrate.CorrectStep(
367 self.cols['CalibData'],
368 self.header['NPIX'],
369 self.header['NROI'],
370 previous_start_cells,
371 self.cols['StartCellData'],
372 3)
373 self.list_of_previous_start_cells.append(self.cols['StartCellData'])
374 if len(self.list_of_previous_start_cells) > 5:
375 self.list_of_previous_start_cells.pop(0)
376
377 for ch in range(self.header['NPIX']):
378 self.drs_calibrate.RemoveSpikes3(self.cols['CalibData2D'][ch], self.header['NROI'])
379
380 return self
381
382
383if __name__ == '__main__':
384 """ Example """
385 print "Example for calibrated raw-file"
386 f = RawData(sys.argv[1], sys.argv[2])
387
388 print "number of events:", len(f)
389 print "date of observation:", f.header['DATE']
390 print "The files has these cols:", f.cols.keys()
391
392 for counter, row in enumerate(f):
393 print "Event Id:", row.cols['EventNum']
394 print "shape of column 'StartCellData'", row.cols['StartCellData'].shape
395 print "dtype of column 'Data'", row.cols['StartCellData'].dtype
396 if counter > 3:
397 break
398 # get next row
399 f.next()
400 print "Event Id:", f.cols['EventNum']
401 # get another row
402 f.next(10)
403 print "Event Id:", f.cols['EventNum']
404 # Go back again
405 f.next(3)
406 print "Event Id:", f.cols['EventNum']
Note: See TracBrowser for help on using the repository browser.