| 1 | #!/usr/bin/python -itt
|
|---|
| 2 | import numpy as np
|
|---|
| 3 | import pprint
|
|---|
| 4 | import ctypes
|
|---|
| 5 |
|
|---|
| 6 | from ROOT import gSystem
|
|---|
| 7 | gSystem.Load('pyfits_h.so')
|
|---|
| 8 | from ROOT import *
|
|---|
| 9 |
|
|---|
| 10 | class FactFits( fits ):
|
|---|
| 11 | """ -Fact Fits File-
|
|---|
| 12 | A Python wrapper for the fits-class implemented in pyfits.h
|
|---|
| 13 | provides easy access to the fits file meta data.
|
|---|
| 14 | * dictionary of file metadata - self.meta
|
|---|
| 15 | * dict of table metadata - self.columns
|
|---|
| 16 | * variable table column access, thus possibly increased speed while looping
|
|---|
| 17 | """
|
|---|
| 18 | def __init__(self, path):
|
|---|
| 19 | """ creates meta and columns dictionaries
|
|---|
| 20 | """
|
|---|
| 21 | self.path = path
|
|---|
| 22 | try:
|
|---|
| 23 | fits.__init__(self,path)
|
|---|
| 24 | except IOError:
|
|---|
| 25 | print 'problem accessing data file: ', data_file_name
|
|---|
| 26 | raise # stop ! no data
|
|---|
| 27 |
|
|---|
| 28 | self.meta = self._make_meta_dict()
|
|---|
| 29 | self.columns = self._make_columns_dict()
|
|---|
| 30 |
|
|---|
| 31 | self.treat_meta_dict()
|
|---|
| 32 |
|
|---|
| 33 |
|
|---|
| 34 | # list of columns, which are already registered
|
|---|
| 35 | # see method register()
|
|---|
| 36 | self._registered_cols = []
|
|---|
| 37 | # dict of column data, this is used, in order to be able to remove
|
|---|
| 38 | # the ctypes of
|
|---|
| 39 | self._table_cols = {}
|
|---|
| 40 |
|
|---|
| 41 | # I need to count the rows, since the normal loop mechanism seems not to work.
|
|---|
| 42 | self._current_row = 0
|
|---|
| 43 |
|
|---|
| 44 | self.stacked_cols = {}
|
|---|
| 45 |
|
|---|
| 46 | def _make_meta_dict(self):
|
|---|
| 47 | """ This method retrieves meta information about the fits file and
|
|---|
| 48 | stores this information in a dict
|
|---|
| 49 | return: dict
|
|---|
| 50 | key: string - all capital letters
|
|---|
| 51 | value: tuple( numerical value, string comment)
|
|---|
| 52 | """
|
|---|
| 53 | # intermediate variables for file metadata dict generation
|
|---|
| 54 | keys=self.GetPy_KeyKeys()
|
|---|
| 55 | values=self.GetPy_KeyValues()
|
|---|
| 56 | comments=self.GetPy_KeyComments()
|
|---|
| 57 | types=self.GetPy_KeyTypes()
|
|---|
| 58 |
|
|---|
| 59 | if len(keys) != len(values):
|
|---|
| 60 | raise TypeError('len(keys)',len(keys),' != len(values)', len(values))
|
|---|
| 61 | if len(keys) != len(types):
|
|---|
| 62 | raise TypeError('len(keys)',len(keys),' != len(types)', len(types))
|
|---|
| 63 | if len(keys) != len(comments):
|
|---|
| 64 | raise TypeError('len(keys)',len(keys),' != len(comments)', len(comments))
|
|---|
| 65 |
|
|---|
| 66 | meta_dict = {}
|
|---|
| 67 | for i in range(len(keys)):
|
|---|
| 68 | type = types[i]
|
|---|
| 69 | if type == 'I':
|
|---|
| 70 | value = int(values[i])
|
|---|
| 71 | elif type == 'F':
|
|---|
| 72 | value = float(values[i])
|
|---|
| 73 | elif type == 'B':
|
|---|
| 74 | if values[i] == 'T':
|
|---|
| 75 | value = True
|
|---|
| 76 | elif values[i] == 'F':
|
|---|
| 77 | value = False
|
|---|
| 78 | else:
|
|---|
| 79 | raise TypeError("meta-type is 'B', but meta-value is neither 'T' nor 'F'. meta-value:",values[i])
|
|---|
| 80 | elif type == 'T':
|
|---|
| 81 | value = values[i]
|
|---|
| 82 | else:
|
|---|
| 83 | raise TypeError("unknown meta-type: known meta types are: I,F,B and T. meta-type:",type)
|
|---|
| 84 | meta_dict[keys[i]]=(value, comments[i])
|
|---|
| 85 | return meta_dict
|
|---|
| 86 |
|
|---|
| 87 |
|
|---|
| 88 | def _make_columns_dict(self):
|
|---|
| 89 | """ This method retrieves information about the columns
|
|---|
| 90 | stored inside the fits files internal binary table.
|
|---|
| 91 | returns: dict
|
|---|
| 92 | key: string column name -- all capital letters
|
|---|
| 93 | values: tuple(
|
|---|
| 94 | number of elements in table field - integer
|
|---|
| 95 | size of element in bytes -- this is not really interesting for any user
|
|---|
| 96 | might be ommited in future versions
|
|---|
| 97 | type - a single character code -- should be translated into
|
|---|
| 98 | a comrehensible word
|
|---|
| 99 | unit - string like 'mV' or 'ADC count'
|
|---|
| 100 | """
|
|---|
| 101 | # intermediate variables for file table-metadata dict generation
|
|---|
| 102 | keys=self.GetPy_ColumnKeys()
|
|---|
| 103 | #offsets=self.GetPy_ColumnOffsets() #not needed on python level...
|
|---|
| 104 | nums=self.GetPy_ColumnNums()
|
|---|
| 105 | sizes=self.GetPy_ColumnSizes()
|
|---|
| 106 | types=self.GetPy_ColumnTypes()
|
|---|
| 107 | units=self.GetPy_ColumnUnits()
|
|---|
| 108 |
|
|---|
| 109 | # zip the values
|
|---|
| 110 | values = zip(nums,sizes,types,units)
|
|---|
| 111 | # create the columns dictionary
|
|---|
| 112 | columns = dict(zip(keys ,values))
|
|---|
| 113 | return columns
|
|---|
| 114 |
|
|---|
| 115 | def stack(self, on=True):
|
|---|
| 116 | self.next()
|
|---|
| 117 | for col in self._registered_cols:
|
|---|
| 118 | if isinstance( self.dict[col], type(np.array('')) ):
|
|---|
| 119 | self.stacked_cols[col] = self.dict[col]
|
|---|
| 120 | else:
|
|---|
| 121 | # elif isinstance(self.dict[col], ctypes._SimpleCData):
|
|---|
| 122 | self.stacked_cols[col] = np.array(self.dict[col])
|
|---|
| 123 | # else:
|
|---|
| 124 | # raise TypeError("I don't know how to stack "+col+". It is of type: "+str(type(self.dict[col])))
|
|---|
| 125 |
|
|---|
| 126 | def register(self, input_str):
|
|---|
| 127 | columns = self.columns
|
|---|
| 128 | if input_str.lower() == 'all':
|
|---|
| 129 | for col in columns:
|
|---|
| 130 | self._register(col)
|
|---|
| 131 | else:
|
|---|
| 132 | #check if colname is in columns:
|
|---|
| 133 | if input_str not in columns:
|
|---|
| 134 | error_msg = 'colname:'+ input_str +' is not a column in the binary table.\n'
|
|---|
| 135 | error_msg+= 'possible colnames are\n'
|
|---|
| 136 | for key in columns:
|
|---|
| 137 | error_msg += key+'\n'
|
|---|
| 138 | raise KeyError(error_msg)
|
|---|
| 139 | else:
|
|---|
| 140 | self._register(input_str)
|
|---|
| 141 |
|
|---|
| 142 | # 'private' method, do not use
|
|---|
| 143 | def _register( self, colname):
|
|---|
| 144 | columns = self.columns
|
|---|
| 145 | local = None
|
|---|
| 146 |
|
|---|
| 147 | number_of_elements = int(columns[colname][0])
|
|---|
| 148 | size_of_elements_in_bytes = int(columns[colname][1])
|
|---|
| 149 | ctypecode_of_elements = columns[colname][2]
|
|---|
| 150 | physical_unit_of_elements = columns[colname][3]
|
|---|
| 151 |
|
|---|
| 152 | # snippet from the C++ source code, or header file to be precise:
|
|---|
| 153 | #case 'L': gLog << "bool(8)"; break;
|
|---|
| 154 | #case 'B': gLog << "byte(8)"; break;
|
|---|
| 155 | #case 'I': gLog << "short(16)"; break;
|
|---|
| 156 | #case 'J': gLog << "int(32)"; break;
|
|---|
| 157 | #case 'K': gLog << "int(64)"; break;
|
|---|
| 158 | #case 'E': gLog << "float(32)"; break;
|
|---|
| 159 | #case 'D': gLog << "double(64)"; break;
|
|---|
| 160 |
|
|---|
| 161 |
|
|---|
| 162 |
|
|---|
| 163 | # the fields inside the columns can either contain single numbers,
|
|---|
| 164 | # or whole arrays of numbers as well.
|
|---|
| 165 | # we treat single elements differently...
|
|---|
| 166 | if number_of_elements == 1:
|
|---|
| 167 | # allocate some memory for a single number according to its type
|
|---|
| 168 | if ctypecode_of_elements == 'J': # J is for a 4byte int, i.e. an unsigned long
|
|---|
| 169 | local = ctypes.c_ulong()
|
|---|
| 170 | un_c_type = long
|
|---|
| 171 | elif ctypecode_of_elements == 'I': # I is for a 2byte int, i.e. an unsinged int
|
|---|
| 172 | local = ctypes.c_ushort()
|
|---|
| 173 | un_c_type = int
|
|---|
| 174 | elif ctypecode_of_elements == 'B': # B is for a byte
|
|---|
| 175 | local = ctypes.c_ubyte()
|
|---|
| 176 | un_c_type = int
|
|---|
| 177 | elif ctypecode_of_elements == 'D':
|
|---|
| 178 | local = ctypes.c_double()
|
|---|
| 179 | un_c_type = float
|
|---|
| 180 | elif ctypecode_of_elements == 'E':
|
|---|
| 181 | local = ctypes.c_float()
|
|---|
| 182 | un_c_type = float
|
|---|
| 183 | elif ctypecode_of_elements == 'A':
|
|---|
| 184 | local = ctypes.c_uchar()
|
|---|
| 185 | un_c_type = chr
|
|---|
| 186 | elif ctypecode_of_elements == 'K':
|
|---|
| 187 | local = ctypes.c_ulonglong()
|
|---|
| 188 | un_c_type = long
|
|---|
| 189 | else:
|
|---|
| 190 | raise TypeError('unknown ctypecode_of_elements:',ctypecode_of_elements)
|
|---|
| 191 | else:
|
|---|
| 192 | if ctypecode_of_elements == 'B': # B is for a byte
|
|---|
| 193 | nptype = np.int8
|
|---|
| 194 | elif ctypecode_of_elements == 'A': # A is for a char .. but I don't know how to handle it
|
|---|
| 195 | nptype = np.int8
|
|---|
| 196 | elif ctypecode_of_elements == 'I': # I is for a 2byte int
|
|---|
| 197 | nptype = np.int16
|
|---|
| 198 | elif ctypecode_of_elements == 'J': # J is for a 4byte int
|
|---|
| 199 | nptype = np.int32
|
|---|
| 200 | elif ctypecode_of_elements == 'K': # B is for a byte
|
|---|
| 201 | nptype = np.int64
|
|---|
| 202 | elif ctypecode_of_elements == 'E': # B is for a byte
|
|---|
| 203 | nptype = np.float32
|
|---|
| 204 | elif ctypecode_of_elements == 'D': # B is for a byte
|
|---|
| 205 | nptype = np.float64
|
|---|
| 206 | else:
|
|---|
| 207 | raise TypeError('unknown ctypecode_of_elements:',ctypecode_of_elements)
|
|---|
| 208 | local = np.zeros( number_of_elements, nptype)
|
|---|
| 209 |
|
|---|
| 210 | # Set the Pointer Address
|
|---|
| 211 | self.SetPtrAddress(colname, local)
|
|---|
| 212 | self._table_cols[colname] = local
|
|---|
| 213 | if number_of_elements > 1:
|
|---|
| 214 | self.__dict__[colname] = local
|
|---|
| 215 | self.dict[colname] = local
|
|---|
| 216 | else:
|
|---|
| 217 | # remove any traces of ctypes:
|
|---|
| 218 | self.__dict__[colname] = local.value
|
|---|
| 219 | self.dict[colname] = local.value
|
|---|
| 220 | self._registered_cols.append(colname)
|
|---|
| 221 |
|
|---|
| 222 |
|
|---|
| 223 | def treat_meta_dict(self):
|
|---|
| 224 | """make 'interesting' meta information available like normal members.
|
|---|
| 225 | non interesting are:
|
|---|
| 226 | TFORM, TUNIT, and TTYPE
|
|---|
| 227 | since these are available via the columns dict.
|
|---|
| 228 | """
|
|---|
| 229 |
|
|---|
| 230 | self.number_of_rows = self.meta['NAXIS2'][0]
|
|---|
| 231 | self.number_of_columns = self.meta['TFIELDS'][0]
|
|---|
| 232 |
|
|---|
| 233 | # there are some information in the meta dict, which are alsways there:
|
|---|
| 234 | # there are regarded as not interesting:
|
|---|
| 235 | uninteresting_meta = {}
|
|---|
| 236 | uninteresting_meta['arraylike'] = {}
|
|---|
| 237 | uninteresting = ['NAXIS', 'NAXIS1', 'NAXIS2',
|
|---|
| 238 | 'TFIELDS',
|
|---|
| 239 | 'XTENSION','EXTNAME','EXTREL',
|
|---|
| 240 | 'BITPIX', 'PCOUNT', 'GCOUNT',
|
|---|
| 241 | 'ORIGIN',
|
|---|
| 242 | 'PACKAGE', 'COMPILED', 'CREATOR',
|
|---|
| 243 | 'TELESCOP','TIMESYS','TIMEUNIT','VERSION']
|
|---|
| 244 | for key in uninteresting:
|
|---|
| 245 | if key in self.meta:
|
|---|
| 246 | uninteresting_meta[key]=self.meta[key]
|
|---|
| 247 | del self.meta[key]
|
|---|
| 248 |
|
|---|
| 249 | # the table meta data contains
|
|---|
| 250 |
|
|---|
| 251 |
|
|---|
| 252 | # shortcut to access the meta dict. But this needs to
|
|---|
| 253 | # be cleaned up quickly!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
|---|
| 254 | meta = self.meta
|
|---|
| 255 |
|
|---|
| 256 | # loop over keys:
|
|---|
| 257 | # * try to find array-like keys
|
|---|
| 258 | arraylike = {}
|
|---|
| 259 | singlelike = []
|
|---|
| 260 | for key in self.meta:
|
|---|
| 261 | stripped = key.rstrip('1234567890')
|
|---|
| 262 | if stripped == key:
|
|---|
| 263 | singlelike.append(key)
|
|---|
| 264 | else:
|
|---|
| 265 | if stripped not in arraylike:
|
|---|
| 266 | arraylike[stripped] = 0
|
|---|
| 267 | else:
|
|---|
| 268 | arraylike[stripped] += 1
|
|---|
| 269 | newmeta = {}
|
|---|
| 270 | for key in singlelike:
|
|---|
| 271 | newmeta[key.lower()] = meta[key]
|
|---|
| 272 | for key in arraylike:
|
|---|
| 273 | uninteresting_meta['arraylike'][key.lower()] = []
|
|---|
| 274 | for i in range(arraylike[key]+1):
|
|---|
| 275 | if key+str(i) in meta:
|
|---|
| 276 | uninteresting_meta['arraylike'][key.lower()].append(meta[key+str(i)])
|
|---|
| 277 | self.ui_meta = uninteresting_meta
|
|---|
| 278 | # make newmeta self
|
|---|
| 279 | for key in newmeta:
|
|---|
| 280 | self.__dict__[key]=newmeta[key]
|
|---|
| 281 |
|
|---|
| 282 | dict = self.__dict__.copy()
|
|---|
| 283 | del dict['meta']
|
|---|
| 284 | del dict['ui_meta']
|
|---|
| 285 | self.dict = dict
|
|---|
| 286 |
|
|---|
| 287 | def __iter__(self):
|
|---|
| 288 | """ iterator """
|
|---|
| 289 | return self
|
|---|
| 290 |
|
|---|
| 291 | def next(self):
|
|---|
| 292 | """ used by __iter__ """
|
|---|
| 293 | # Here one might check, if looping makes any sense, and if not
|
|---|
| 294 | # one could stop looping or so...
|
|---|
| 295 | # like this:
|
|---|
| 296 | #
|
|---|
| 297 | # if len(self._registered_cols) == 0:
|
|---|
| 298 | # print 'warning: looping without any registered columns'
|
|---|
| 299 | if self._current_row < self.number_of_rows:
|
|---|
| 300 | if self.GetNextRow() == False:
|
|---|
| 301 | raise StopIteration
|
|---|
| 302 | for col in self._registered_cols:
|
|---|
| 303 | if isinstance(self._table_cols[col], ctypes._SimpleCData):
|
|---|
| 304 | self.__dict__[col] = self._table_cols[col].value
|
|---|
| 305 | self.dict[col] = self._table_cols[col].value
|
|---|
| 306 |
|
|---|
| 307 | for col in self.stacked_cols:
|
|---|
| 308 | if isinstance(self.dict[col], type(np.array(''))):
|
|---|
| 309 | self.stacked_cols[col] = np.vstack( (self.stacked_cols[col],self.dict[col]) )
|
|---|
| 310 | else:
|
|---|
| 311 | self.stacked_cols[col] = np.vstack( (self.stacked_cols[col],np.array(self.dict[col])) )
|
|---|
| 312 | self._current_row += 1
|
|---|
| 313 | else:
|
|---|
| 314 | raise StopIteration
|
|---|
| 315 | return self
|
|---|
| 316 |
|
|---|
| 317 | def show(self):
|
|---|
| 318 | pprint.pprint(self.dict)
|
|---|
| 319 |
|
|---|
| 320 | if __name__ == '__main__':
|
|---|
| 321 | import sys
|
|---|
| 322 | if len(sys.argv) == 1:
|
|---|
| 323 | print 'usage:', sys.argv[0], 'fits-file-name'
|
|---|
| 324 |
|
|---|
| 325 | file = FactFits(sys.argv[1])
|
|---|
| 326 | print '-'*70
|
|---|
| 327 | print "opened :", sys.argv[1], " as 'file'"
|
|---|
| 328 | print
|
|---|
| 329 | print '-'*70
|
|---|
| 330 | print 'type file.show() to look at its contents'
|
|---|
| 331 | print "type file.register( columnname ) or file.register('all') in order to register columns"
|
|---|
| 332 | print
|
|---|
| 333 | print " due column-registration you declare, that you would like to retrieve the contents of one of the columns"
|
|---|
| 334 | print " after column-registration, the 'file' has new member variables, they are named like the columns"
|
|---|
| 335 | print " PLEASE NOTE: immediatly after registration, the members exist, but they are empty."
|
|---|
| 336 | print " the values are assigned only, when you call file.next() or when you loop over the 'file'"
|
|---|
| 337 | print
|
|---|
| 338 | print "in order to loop over it, just go like this:"
|
|---|
| 339 | print "for row in file:"
|
|---|
| 340 | print " print row.columnname_one, row.columnname_two"
|
|---|
| 341 | print
|
|---|
| 342 | print ""
|
|---|
| 343 | print '-'*70
|
|---|
| 344 | |
|---|