readMilleBinary.py 6.2 KB
Newer Older
1
#!/usr/bin32/python
2
3
4
5

## \file
# Read millepede binary file and print records
#
Claus Kleinwort's avatar
Claus Kleinwort committed
6
# \author Claus Kleinwort, DESY, 2009 (Claus.Kleinwort@desy.de)
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
#
#  \copyright
#  Copyright (c) 2009 - 2015 Deutsches Elektronen-Synchroton,
#  Member of the Helmholtz Association, (DESY), HAMBURG, GERMANY \n\n
#  This library is free software; you can redistribute it and/or modify
#  it under the terms of the GNU Library General Public License as
#  published by the Free Software Foundation; either version 2 of the
#  License, or (at your option) any later version. \n\n
#  This library is distributed in the hope that it will be useful,
#  but WITHOUT ANY WARRANTY; without even the implied warranty of
#  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#  GNU Library General Public License for more details. \n\n
#  You should have received a copy of the GNU Library General Public
#  License along with this program (see the file COPYING.LIB for more
#  details); if not, write to the Free Software Foundation, Inc.,
#  675 Mass Ave, Cambridge, MA 02139, USA.
#
24
25
# Hardcoded defaults can be replaced by command line arguments for
#    -  Name of binary file
26
#    -  Number of records to print (-1: all; <-1: all, record headers only)
27
#    -  Number of records to skip (optional)
Claus Kleinwort's avatar
Claus Kleinwort committed
28
#    -  Mininum value to print derivative
29
30
31
32
33
34
35
36
37
38
39
40
41
#
# Description of the output from readMilleBinary.py
#    -  Records (tracks) start with \c '===' followed by record number and length 
#       (<0 for binary files containing doubles)
#    -  Measurements: A measurement with global derivatives is called a 'global measurement', 
#       otherwise 'local measurement'. Usually the real measurements from the detectors are 'global'
#       ones and virtual measurements e.g. to describe multiple scattering are 'local'.
#    -  'Global' measurements start with \c '-g-' followed by measurement number, first global label,
#       number of local and global derivatives, measurement value and error. The next lines contain 
#       local and global labels (array('i')) and derivatives (array('f') or array('d')).
#    -  'Local' measurements start with \c '-l-' followed by measurement number, first local label, 
#       number of local and global derivatives, measurement value and error. The next lines contain
#       local labels (array('i')) and derivatives (array('f') or array('d')).
42
#
43
44
# Tested with SL4, SL5, SL6

Claus Kleinwort's avatar
Claus Kleinwort committed
45
import array, sys
46
47
48
49
50
51
52
53
54
55

# ############### read millepede binary file #################
#  
## Binary file type (C or Fortran)
Cfiles = 1  # Cfiles
#Cfiles = 0 # Fortran files
# 
## Integer format 
intfmt = 'i'  # SL5, gcc-4
#intfmt = 'l' # SL4, gcc-3
56
#
57
## Binary file name
Claus Kleinwort's avatar
Claus Kleinwort committed
58
fname = "milleBinaryISN.dat"
59
#
60
## number of records (tracks) to show
61
mrec = 10
62
## number of records (track) to skip before 
63
skiprec = 0
Claus Kleinwort's avatar
Claus Kleinwort committed
64
## minimum value to print derivatives
65
minval = None  # allows for NaNs
66
#
67
# ## C. Kleinwort - DESY ########################
68

69
# ## use command line arguments ?
Claus Kleinwort's avatar
Claus Kleinwort committed
70
71
72
narg = len(sys.argv)
if narg > 1:
  if narg < 3:
73
    print " usage: readMilleBinary.py <file name> <number of records> [<number of records to skip> <minimum value to print derivative>]"
Claus Kleinwort's avatar
Claus Kleinwort committed
74
75
76
77
78
79
    sys.exit(2)
  else:
    fname = sys.argv[1]
    mrec = int(sys.argv[2])
    if narg > 3:
      skiprec = int(sys.argv[3])
Claus Kleinwort's avatar
Claus Kleinwort committed
80
81
    if narg > 4:
      minval = float(sys.argv[4])
82

Claus Kleinwort's avatar
Claus Kleinwort committed
83
#print " input ", fname, mrec, skiprec
84

85
f = open(fname, "rb")
Claus Kleinwort's avatar
Claus Kleinwort committed
86

87
nrec = 0
88
try:
Claus Kleinwort's avatar
Claus Kleinwort committed
89
90
    while (nrec < mrec + skiprec) or (mrec < 0):
# read 1 record
91
92
93
94
95
        nr = 0
        if (Cfiles == 0):
            lenf = array.array(intfmt)
            lenf.fromfile(f, 2)

96
97
98
99
        length = array.array(intfmt)
        length.fromfile(f, 1)
        nr = abs(length[0] / 2)
        nrec += 1
100

101
        if length[0] > 0:
102
            glder = array.array('f')
Claus Kleinwort's avatar
Claus Kleinwort committed
103
        else:
104
            glder = array.array('d')
105
        glder.fromfile(f, nr)
106

107
108
        inder = array.array(intfmt)
        inder.fromfile(f, nr)
109
110
111
112

        if (Cfiles == 0):
            lenf = array.array(intfmt)
            lenf.fromfile(f, 2)
113

114
        if (nrec <= skiprec):  # must be after last fromfile
115
116
            continue

117
        print " === NR ", nrec, length[0] / 2
118

119
120
121
122
        # no details, only header
        if (mrec < -1):
            continue

123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
        i = 0
        nh = 0
        ja = 0
        jb = 0
        jsp = 0
        nsp = 0
        while (i < (nr - 1)):
            i += 1
            while (i < nr) and (inder[i] != 0): i += 1
            ja = i
            i += 1
            while (i < nr) and (inder[i] != 0): i += 1
            jb = i
            i += 1
            while (i < nr) and (inder[i] != 0): i += 1
            i -= 1
139
# special data ?
140
            if (ja + 1 == jb) and (glder[jb] < 0.):
141
142
143
144
145
146
                jsp = jb
                nsp = int(-glder[jb])
                i += nsp
                print ' ### spec. ', nsp, inder[jsp + 1:i + 1], glder[jsp + 1:i + 1]
                continue
            nh += 1
147
            if (jb < i):
148
# measurement with global derivatives
149
                print ' -g- meas. ', nh, inder[jb + 1], jb - ja - 1, i - jb, glder[ja], glder[jb]
150
151
            else:
# measurement without global derivatives
152
                print ' -l- meas. ', nh, inder[ja + 1], jb - ja - 1, i - jb, glder[ja], glder[jb]
153
            if (ja + 1 < jb):
154
155
156
157
158
159
160
161
162
163
164
                lab = []
                val = []
                for k in range(ja + 1, jb):
                    if minval is None:
                        lab.append(inder[k])
                        val.append(glder[k])
                    elif abs(glder[k]) >= minval:
                        lab.append(inder[k])
                        val.append(glder[k])
                print " local  ", lab
                print " local  ", val
165
            if (jb + 1 < i + 1):
166
167
168
169
170
171
172
173
174
175
176
177
178
                lab = []
                val = []
                for k in range(jb + 1, i + 1):
                    if minval is None:
                        lab.append(inder[k])
                        val.append(glder[k]) 
                    elif abs(glder[k]) >= minval:
                        lab.append(inder[k])
                        val.append(glder[k])
                print " global ", lab
                print " global ", val


179
except EOFError:
180
181
    print
    if (nr > 0):
Claus Kleinwort's avatar
Claus Kleinwort committed
182
        print " >>> error: end of file before end of record", nrec
183
184
    else:
        print " end of file after", nrec, "records"
Claus Kleinwort's avatar
Claus Kleinwort committed
185

186
187

f.close()