Version 2.4
25/01/2005
Bug fixes/changes
A reconstruction package for MPGD TPC data. Designed for flexibility, not for speed. With only a few adaptions you should be able to reconstruct and fit tracks. The data can be multiplexed, pads in random order.
Requirements:
F or f95 compiler
HBOOK and MINUIT, i.e f77 CERN libraries
A TAR file for download. It will create
two directories:
mycode containes the library, the source code, and a make file.
work contains user specific files: input files, MYanalyzer.f95 and
a make file which creates TPCanalyzer.run.
In the makefiles change the variables for the directories to your settings.
The determination of the amplitude depends on your electronic. Therefore you
will have to check/make changes in TPCcalcAmplitude.f95.
To read .slcio files you have to install LCIO.
There are three input files you will always
need: TPClayout.txt, TPCfiles.txt and TPCinput.txt.
Optionally there are additional input files to provide information
for each group on falltime, risetime, gain and pedestals as well as
bias correction for rows.
| Supported data formats: the suffix is recognised to read the data in the correct format. | |
| suffix | data format |
| .mid | MIDAS data |
| .mc | JTPC Monte Carlo |
| .dd | Dense data |
| .slcio | LCIO |
| Modules for global variables: | |
| IOunits.f95 | unit numbers for IO files |
| Threshold.f95 | steering parameters, thresholds, ... |
| TPCdata.f95 code |
main data module USE: IOunits, Threshold, Space |
| Modules containing the code: | |
| MYanalyzer.f95 code |
USER program (event loop) USE: IOunits, Threshold, TPCdata, TPCpaw, TPCcalcAmplitude, TPCtrack, Space |
| TPCanalyzer.f95 code |
MAIN program and analysis flow USE: IOunits, Threshold, TPCdata, MYanalyzer, TPCpaw, TPCio, TPCcalcAmplitude |
| TPCcalcAmplitude.f95 code |
calculate signal amplitude from ADC spectrum USE: IOunits, Threshold, TPCdata, TPCpaw, MyMinuit, Space |
| TPCtrack.f95 code |
event reconstruction, tracking USE: IOunits, Threshold, TPCdata, TPCpaw, MyMinuit, Space |
| TPCio.f95 code |
input / output USE: IOunits, Threshold, Space, TPCdata, DataFormats, TPCpaw, MyMinuit |
| DataFormats.f95 code |
knows how to deal with different data format USE: IOunits, TPCdata, Threshold, Space |
| LCIO4TPC.F95 code |
wrapper for LCIO USE: LCIOapi, TPCdata |
| LCIOdummy.f95 code |
replacement for LCIO4TPC.F95 if LCIO is not installed |
| LCIOapi.F95 code |
header file for LCIO routines, parameters |
| TPCpaw.f95 code |
booking and filling of Ntuples and histograms USE: IOunits, Threshold, TPCdata, Space |
| TPCf77.f code |
f77 help routines |
| Self contained code Modules: | |
| MyMinuit.f95 code |
wrapper for MINUIT |
| Space.f95 code |
knows what points, lines, tracks, rectangles, etc. are |
| Input files: | |
| TPCfiles.txt txt |
which MIDAS files will be processed |
| TPCinput.txt txt |
steering cards, thresholds, ... |
| TPClayout.txt | pad layout | TPCpedestal.txt | pedestals for each group |
| TPCftime.txt | fall times for each group |
| TPCgain.txt | gain corrections for each group |
TPCdata.f95
contains the complete data structure and procedures to fill and retrieve signal
amplitudes.
Procedures:
| ReadChannel: | reads the layout from file |
| FillSignal: | fill amplitude and T0 of a pad. The amplitude must be larger than minSignal. T0 is between 1 and nTbin. |
| VoidSignal: | sets the amplitude of a pad to -1. |
| PadGroup: | returns the group of a pad |
| PadRow: | returns the row of a pad |
| PadLocation: | returns the location of a pad (one of 2 places where the location_type of the pad is hard coded). |
| getSignal: | get amplitude and T0 of a pad. |
| getSignalGroup: | get amplitude and T0 of a group. In principle different pads of a group can have different amplitudes. The maximum amplitude and the time of this pad in returned. |
| overflowGroup: | number of timebins with ADC overflow in group i. |
| getSignalRow: | get amplitude and T0 of a row. Sum the amplitudes of groups in this row, T0 is the amplitude weighted average. |
| GmaxSignalRow: | get group with maximum amplitude in a row, return its amplitude and T0. |
| Neighbour: | gets the pad numbers to the left and right. |
| AllocArrays: | allocates the arrays listed below. |
Parameters:
| Layout: | |
| For multiplexing pads are combined in groups which correspond to the readout channels. Pads with the same value for the Y position are combined in rows. Multiplexing across rows is possible. Unless for track fitting rows are not required, i.e. each pad can have its own row number. | |
| nPads,nGroups,nRows,nTbin: | number of pads, groups, rows and time bins. |
| maxPinGroup,maxPinRow,maxGinRow: | maximum number of pads-in-groups, pads-in-rows, groups-in-rows. |
| nPadsinGroup(nGroups): | actual number of pads in group i. |
| nPadsinRow(nRows): | actual number of pads in row i. |
| nGroupsinRow(nRows): | actual number of groups in row i. |
| iPadofGroup(maxPinGroup,nGroups): | list of pads in group i. |
| iPadofRow(maxPinRow,nRows): | list of pads in row i. |
| iGroupofRow(maxGinRow,nRows): | list of groups in row i. |
| NumMaps: | Number of different pad shapes for Grid_type. |
| Event data: | |
| RunNumber, EventNumber, Time | |
| IADC(nTbin,nGroups): | ADC spectrum of nTbin time bins for all groups. |
| ADC(nTbin,nGroups): | Pedestal subtracted ADC spectrum. |
| MCx,MCz0,MCphi,MCtheta: | Monte Carlo Tree information. |
| info(nPads) (PRIVATE): | information for all pads, is of type Pad_type |
| Pad_type (PRIVATE) : | contains group, row, amplitude, T0 and location location is of a type characterising the pad, i.e. Rectangle_type or Grid_type. (One of 2 places in the code where the location type of the pad is specified. Change it here and the rest takes (hopefully) care of itself). |
MYanalyzer.f95
USER program (event loop)
| analyzEvent: | calculate signal amplitudes, correct gain, preselection, find seed track, (demultiplex), (fit all pulses), track fit Z-Y, track fit X-Y, resolution fits, fill analysis Ntuple (or whatever you want) |
| MY_book: | booking of user histograms and ntuples. |
| NF_analysis (PRIVATE): | Example of user Ntuple. |
| HF_test (PRIVATE): | Example of user histogram. |
| DenseData: | Write DenseData. |
TPCanalyzer.f95
MAIN program and analysis flow
| PROGRAM TPCanalyzer: | initialize, loop runfiles, loop events, read data, call analyzEvent (or Calibrate), terminate program. |
TPCcalcAmplitude.f95
Calculate signal amplitude from ADC spectrum
| Calibration: pedestals, gain, falltime: | |
| Calibrate: | no event selection, for all events fit all pulses (with large amplitude), fill calibration histograms. calcPedestal is called or default values are used. |
| correctGain: | divides the amplitudes of all pads by a factor dependent on the group. |
| calcPedestal: | calculate pedestal as average over 5 time bins, the maximum of the first or last time bins is taken. Set default values for gain, rise- and fall time. |
| subPedestal: | fills ADC() = IADC() - pedestal. |
| Several routines calculate the signals for all pads: | |
| calcAmplitudePad: | Independent determination for each pad. Determine TMax as time bin with maximum amplitude, calculate amplitude with routines given below, calculate T0, correct amplitude from TMax to T0. |
| calcAmplitudeRow: | Determination within rows. 1st calculate amplitudes and T0 independently with calcAmplitudePad; an amplitude weighted average T0ROW is determined, for pulses <smallSignal the amplitude is recalculated wrt T0ROW. |
| calcAmplitudeFoil: | T0 is determined from largest pulse in row, amplitude is the average over 30 to 100 timebins, depending on risetime of the pulse. |
| calcAmplitudeFoilNEW: | T0 is determined from largest pulse in row, amplitude is the average over 100 timebins using the same window for all pulses in the row. |
| Several routines calculate the amplitude from ADC spectrum depending on a given T (PRIVATE): | |
| calcSignalFix: | Subtract pedestal and return amplitude = -ADC count at given T; to be used only for testing. |
| calcSignalSum: | amplitude is average in range [T0-riseTime,T0+fallTime]. |
| calcSignalSlideSum: | amplitude is average in range [T0-nTsum/2,T0+nTsum/2], nTsum is argument. If T0 is not given, it will be |
| calcSignalExp: | call amplitudeExp, determine T0 via calcT0, correct amplitude to T0. |
| calcSignalExpT0: | call amplitudeExp, correct amplitude to given T0. |
| amplitudeExp: | Subtract pedestal, integrate ADC spectrum for time bins T+50 to T+350, if the pulse is late the window is decreased down to 100 time bins; T should be TMax. |
| 2 routines calculate the time of the pulse (PRIVATE): | |
| calcTMax: | time where the ADC value is largest. |
| calcT0: | time of the rising pulse where the ADC value is half the amplitude at TMax. |
| Pulse fitting: | |
| fitAllPulses: | CPU time consuming, should be called only for good events. 700 time bins are used (200 before T0) to fit pulseShape to the ADC spectrum. Start values taken from calcAmplitudeXxx, 2 possibilities: calibration run: fit only large pulses (amplitude >20), 5 free parameters: amplitude, T0, rise time, fall time, pedestal. VERY CPU time consuming. Fit results for rise time, fall time and pedestal are filled in histograms. The average will be written in calibration files. normal run: fit pulses with non-zero amplitudes (>2), 2 free parameters: amplitude, T0. T0 is bound within +- 200 bins of the start value. |
| Plotting: books and fills histograms with pulse shapes | |
| plotPulseAll: | all groups. |
| plotPulseEvent: | groups with amplitude > 0. |
| plotPulseRow: | using common time in one row. |
| plotPulseVeto: | for groups given as argument. |
TPCtrack.f95
Event reconstruction, tracking
| eventVeto: | veto events if amplitudes in veto groups (iGroupV) too high, not enough hits in seed rows (iRow1,iRow2), not enough hits in track rows (iRowT). |
| findSeedTrack: | find a seed track to be used as start value for the track fit. There are two sets of rows given (iRow1,iRow2). In each set the pad with the highest amplitude is found. A X-Y coordinate is determined as an amplitude weighted average of this pad and its neighbours within the row-set. The line connecting the coordinates in the two sets is the seed track. |
| deMultiplex: | set the amplitudes of all pads to -1 apart from the one pad in each group that is closest to the seed track. |
| time (Z-Y) track fits: | |
| TimeFit: | perform a Z-Y track fit, 2 parameters: T(Y=0) and a slope, track rows iRowT are used. |
| TimeLine: | equivalent to TimeFit but faster, uses linear regression. |
| TimeFitZ: | same as TimeLine but returns a line (Ztrack). |
| Spatial (X-Y) track fits: Overview. The seedTrack is input and the fitted track is returned on the same parameter. The type can be a line or a XYtrack. In the latter case the radius is fitted too. The likelihood is determined from a uniform line of charge with gaussian profile or PadResponsFunction (LPRF=.TRUE.), expected amplitudes are normalised within each row, negative loglikelihood is minimized with MINUIT. | |
| TrackFitWidth: | perform a X-Y track fit, width is free parameter, additional parameters: X(Y=0), PHI (and radius). |
| TrackFitXY: | perform a X-Y track fit, width is fixed depending on drift distance z, parameters: X(Y=0), PHI (and radius). |
| TrackFitX: | only X(Y=0) is free, the other parameters are fixed to the values of the track given. |
| ROWCHI (PRIVATE): | likelihood for trackfit. Takes all pads in the row into account. |
| ROWCHItrunc (PRIVATE): | likelihood for trackfit. Uses only a given number of pads closest to the track. |
| ROWCHIPRF (PRIVATE): | chisq for trackfit. |
| AERR: | Error on amplitudes for ROWCHIPRF. Needs parameters AERRpar be set properly. |
| PRF: | calculates the integrated charge on a pad from Gaussian charge profile or pad response function, takes bias into account. |
| QRAT (PRIVATE): | ratio of quartics for pad response function. Needs parameters TrackWidth, PRFdel, PRF be set properly. |
TPCio.f95
input / output
| IOinitialize: | open all files, read which data to process (TPCfiles.txt), initialize PAW and MINUIT, call ReadInput and ReadLayout. |
| IOend: | close all files, write histograms and Ntuples to file. |
| OPENData: | open the data file, determine from the extention which data format to use. |
| ReadHeader: | wrapper to read the run header in the right format. |
| ReadEvent: | wrapper to read the event in the right format. |
| readCalibData: | for normal runs reads information from TPCcalib.txt, TPCpedestal.txt, TPCgain.txt, TPCftime.txt. If files are not available defaults are used. |
| ReadInput (PRIVATE): | read information stored in Threshold from TPCinput.txt. |
| ReadLayout (PRIVATE): | read pad layout from TPClayout.txt, allocate arrays. |
DataFormats.f95
knows the different data formats
| Write DenseData: | |
| WriteDenseHeader: | write run header. |
| WriteDenseData: | write event, optionally including additional information (array), e.g. track. |
| WriteEndOfRun: | write end of run mark. |
| Read DenseData: | |
| OPENDense: | open DenseData file. |
| ReadDenseHeader: | read run header. |
| ReadDenseEvent: | read event. |
| JTPC Monte Carlo: | |
| ReadMCHeader: | read run header. |
| ReadMCEvent: | read event. |
| MIDAS: | |
| OPENMIDAS: | open MIDAS or JTPC file. |
| ReadMIDASHeader: | read run header, fills number of time bins and number of groups. |
| ReadMIDASEvent: | read event, fill ADC spectrum for all groups. |
| LCIO: | |
| OPENLCIO: | open LCIO file. |
| ReadLCIOHeader: | read comment from run header, scan first 10 events for information on time bins. |
| ReadLCIOEvent: | read event and MCparticle if available. |
LCIO4TPC.F95
Wrapper for LCIO FORTRAN interface
LCIO parameters are defined in LCIOapi.F95.
| OpenReadLCIO: | open LCIO file for reading. |
| ReadLCIOnTbin: | if LCIO_TPCBIT_RAW bit is set (raw data available) search first 10 events for information on nTbin. |
| ReadLCIOComment: | read comment from run header. |
| CloseReadLCIO: | close LCIO file. |
| ReadLCIOnextEvent: | get next event; read Amplitude/Time or raw data; read MCtrack
if available. will read raw data if available unless bit 1 in iflag is set, otherwise will read Amplitude/Time. |
| OpenWriteLCIO: | open LCIO file for writing; file name = thisrun.slcio. |
| CloseWriteLCIO: | close LCIO file. |
| WriteLCIOHeader: | write comment / date to run header. |
| WriteLCIOEvent: | depending whether bit in iflag is set: bit 0: write raw data (if nTbin > 0) bit 1: write Amplitude / Time . |
LCIOdummy.F95
Replacement for LCIO4TPC.F95 if LCIO is not installed
Same procedures available, but they do nothing.
LCIOapi.F95
Definitions of LCIO functions and parameters (headerfile)
To replace the include files of f77. Not complete yet. For functions with side effects (result returned as formal parameter) wrapper routines (function or subroutine) are defined. Naming convention get -> got.
TPCpaw.f95
Booking and filling of Ntuples and histograms
Try to keep PAW histogram IDs in one place.
| PAWini: PAWend: |
initilize and terminate PAW. |
| PAW_book : | book all histograms and Ntuples. |
| NT_book: | wrapper to book ntuple Ntuple. |
| NF_xxxxx: | fill Ntuple, e.g. NF_event for event display. |
| HF_xxxxx: | fill histogram. |
| PininCal: | book histograms for calibration: pedestals, ... |
| calSignal: | average amplitude, pedestal; fill calibration histograms. |
| calRFtime: | rise- and fall-time; fill calibration histograms. |
TPCf77.f
f77 help routines
| OPENFILE: CLOSEFILE: |
wrapper for OPEN and CLOSE (f77 can't write to a file opened with F!). |
| WRITECH: | write a character string to a file. |
| PAWCOM: | initialize PAW, PAW common block. |
| READWORD, READHWORD, READBYTE, READCH: | |
| wrapper for FGETC for MIDAS data. | |
MyMinuit.f95
Wrapper for MINUIT
| MINUIT: | perform fit. input: subroutine that calculates the chisq (can be in F), arrays containing names, start values, errors and bounds for each fit parameter, print flag. output: fit values, errors, error matrix, chisq, error flag. |
| MINini: MINend: |
initialize and terminate MINUIT. |
| LinearRegression: | linear fit via linear regession. |
| Qfit: | quartic fit via differentiation and linear regession. |
| SPLINE: | spline intrapolation. |
Space.f95
knows what points, lines, tracks, rectangles, etc. are
Contains procedures to fill objects and get parameters back.
| Coordinate_type: | X and Y coordinate. |
| Line_type: | X0 and PHI, can be filled with these parameters or any coordinate and PHI or two coordinates (line connecting two points). Also tangent of a XYtrack at given y. |
| XYTrack_type: | a curved line, to be filled with a line (tangent at Y=0) and the inverse radius. |
| Rectangle_type: | X,Y coordinate for centre and X,Y coordinate for pitch(=width) and length(=height). |
| Grid_type: | a rectangle subdivided in nX times nY small rectangles (recty), a bitmap specifies which recty belongs to the pad. |
| getCentre: | centre (coordinate) of a Pad_type. |
| Distance: | calculate the distance between 2 objects. |
| FillLoc: | fill the location of a Pad_type. |
| AllocBitMap: | allocates the bit map for Grid_type. |
| FillBitMap: | fills the bit map for Grid_type. |
| printGrid: | prints the grid information for Grid_type. |
| intTrackPad: | integral of a track over a pad. |
| operators: + - * == = |
partially defined for coordinates and lines. |
| TPCfiles.txt: | which MIDAS files will be processed. F number: number of MIDAS files to be read, must be given. R nev1 nev2: range of events to be processed, running event number. S nev1 nev2 nev3 ...: sample of events to be processed, MIDAS event number. If neither R nor S are given the full run is processed; range and sample do not work if more than one file shall be processed. *: comment, line ignored. > MIDAS files: a list of MIDAS files should follow this line. |
| TPCinput.txt: | steering cards, thresholds, ... You should know how these parameters are used. There are default values for all of them but they might be useless. You might not use all parameters for your analysis. > keyword: the following line contains the value for keyword. *: comment, line ignored. |
| TPClayout.txt: | pad layout. top: values for nPads,nGroups,nRows,nTbin,maxPinGroup,maxPinRow,maxGinRow For rectangles: list of pad, group, row, centre position along the row (along width), centre-position of the row (along height), width and height for each pad. For Grid_type in addition a grid number is given. The grid information is given at the bottom: Grid number, number of columns and rows in one line followed by the pattern. space and . mean empty, every other sign means 'belongs to pad'. This is also an example for multiplexed layout. Apart from the pattern everything is read without format (i.e. spacing does not matter). But there are no real comment lines, the number of lines skipped is fixed. An example how to write the layout file. The pads, groups and rows must be given from 1 to N, pads must be in order. If your input data has groups (read-out channels) in the range from 0 to N-1, group 0 will be called N when the data are read in. Please change your layout accordingly. |
| TPCpedestal.txt: | pedestals for each group. One comment line followed by a line for each group which gives the value. |
| TPCftime.txt: | fall times for each group; format as for pedestals. |
| TPCrtime.txt: | rise times for each group; format as for pedestals. |
| TPCgain.txt: | relative gain for each group; format as for pedestals. |
| TPCbias.txt: | bias for rows. 1st line: comment 2nd line: numbers of rows nR for which a bias is given, number of x-values nX (should be the same for all rows). 3rd line: list of nR row numbers, the first number is a dummy. nX lines with values for: X (position within row), nR bias values. |