MyMarlinTPC
170316
|
Simultaneous hit and track finding (pad row based). More...
#include <RowBasedPadPulseRoadSearchProcessor.h>
Public Member Functions | |
virtual Processor * | newProcessor () |
RowBasedPadPulseRoadSearchProcessor () | |
Construct processor. More... | |
virtual void | init () |
Initialize processor. More... | |
virtual void | processRunHeader (EVENT::LCRunHeader *run) |
virtual void | processEvent (EVENT::LCEvent *evt) |
Process event. More... | |
virtual void | check (EVENT::LCEvent *evt) |
virtual void | end () |
Protected Attributes | |
std::string | _inputTrackerPulsesColName |
Name of the TPc pulse input collection. More... | |
std::string | _outputTrackerHitsColName |
Name of the TPC hit output collection. More... | |
std::string | _outputTrackColName |
Name of the track output collection. More... | |
double | _bfieldScaleFactor |
scale factor for magnetic field (default: 1.0) More... | |
double | _driftVelocity |
drift velocity More... | |
int | _maxColDiff |
maximal column difference to pulse in neighbourhood More... | |
double | _dtMax |
maximum drift time difference in neighboorhood/hit More... | |
int | _maxSeedsInRow |
maximum number of seeding pulses in row More... | |
int | _minRowDiff |
minimum row difference for pulse pair seeding road More... | |
double | _maxDxyRoad |
maximum interpolation deviation in XY (road search) More... | |
double | _maxDzRoad |
maximum interpolation deviation in Z (road search) More... | |
int | _minRows |
minimum number of rows for segment seed More... | |
double | _rowDensityCut |
minimum row 'density' for segment seed More... | |
double | _maxVarXY |
maximum (average) XY variance (unweighted Chi2/ndf) for segment seed More... | |
double | _maxOverlapFraction |
maximum overlap fraction of segment seeds More... | |
int | _maxColDiffHit |
maximum distance of pad to seed in hit (if > 0) More... | |
int | _minNumPulses |
minimum number of pulses in hit More... | |
int | _maxColGap |
maximum column gap in hit More... | |
double | _offsetCol |
measurement variance: sigma0^2 in column direction More... | |
double | _slopeCol |
measurement variance: diffusion term in column direction More... | |
double | _offsetDrift |
measurement variance: sigma0^2 in drift direction More... | |
double | _slopeDrift |
measurement variance: diffusion term in drift direction More... | |
double | _offsetRow |
measurement variance: sigma0^2 in row direction (use rowHeight^2/12. if < 0.) More... | |
double | _slopeRow |
measurement variance: diffusion term in row direction More... | |
double | _chargeConversionFactor |
charge conversion factor from ADC values to primary electrons More... | |
double | _segCut |
Chi2/Ndf cut for segment matching (Ndf=4 for B off, 5 for B on) More... | |
bool | _refAtPCA |
Use Pca as reference point (else 1. hit) More... | |
bool | _skipRoadSearch |
skip road search, hit finding only More... | |
bool | _skipMultiplePulses |
skip multiple pulse candidates as road search seeds More... | |
Private Member Functions | |
bool | areNeighbourModules (gear::TPCModule *, gear::TPCModule *, bool) |
Check if modules are neighbours. More... | |
Private Attributes | |
double | _Bzc |
magnetic field strength (Bz*c) More... | |
std::map< int, std::vector< int > > | _modNeighbours |
neighbour modules More... | |
double | _Xcenter |
TPC center, X coordinate. More... | |
double | _Ycenter |
TPC center, Y coordinate. More... | |
std::map< int, int > | _moduleRowOffset |
row offsets per module layer More... | |
Simultaneous hit and track finding (pad row based).
Input are the (position and measurement direction) of the TPC pad pulses, a parametrization of the hit measurement errors and module parameters (offset and extend) from GEAR.
In each module (with pulses):
Segments from each module are combined with compatible segments from the corresponding neighbouring modules to build track candidates (using only unique matches). As reference point for the track candidates the first hit is used (optional: PCA).
More details have been presented at the DESY FLC/TPC meeting on the 9th of February 2017.
The hits on the track candidates are ordered by row number, no multiple hits in a row are allowed. All matching criteria use chi2 cuts based on the hit measurement errors. The calculation of the measurement (pad) direction in XY uses "module.getLocalPadLayout().getCoordinateType()" to select polar or cartesian pad row geometry.
The default parameters value in '()' below correspond to B=1.0 T data, alternative settings for B=0 are given in '[]'.
"InputTrackerPulses":string | The name of the input collection of TPC (pad) pulses ("TPCPulses") |
"OutputTrackerHits":string | The name of the output collection with the found hits ("TPCHits") |
"OutputTracks":string | The name of the output collection with the found tracks ("TPCTracks") |
"BFieldScaleFactor":double | Optional parameter, scales magnetic field (map), use 1.0 (default) for field ON or 0.0 for field OFF (1.) [0.] |
"DriftVelocity":double | Optional parameter, drift velocity (74.) |
"MaxColDiffNeighbourhood":int | Optional parameter, maximal column difference to pulse in neighbourhood (1) [2] |
"MaxTimeDiff":double | Optional parameter, maximum drift time difference in neighboorhood/hit (80.) |
"MaxSeedsInRow":int | Optional parameter, maximum number of seeding pulses in row (10) |
"MinRowDiff":int | Optional parameter, minimum row difference for pulse pair seeding road (5) |
"MaxDxyRoad":double | Optional parameter, maximum interpolation deviation in XY (road search) (2.) [4.] |
"MaxDzRoad":double | Optional parameter, maximum interpolation deviation in Z (road search) (3.) |
"MinRowNum":int | Optional parameter, minimum number of rows for segment seed (6) |
"MinRowDensity":double | Optional parameter, minimum row 'density' for segment seed (0.8) |
"MaxVarXY":double | Optional parameter, maximum (average) XY variance (unweighted Chi2/ndf) for segment seed (0.25) [2.] |
"MaxOverlapFraction":double | Optional parameter, maximum overlap fraction of segment seeds (0.5) |
"MaxColDiffHit":int | Optional parameter, maximum distance of pad to seed in hit (if > 0) (0) |
"MinNumPulses":double | Optional parameter, minimum number of pulses in hit (2) |
"MaxColGap":double | Optional parameter, maximum column gap in hit (1) |
"OffsetCol":double | Optional parameter, hit variance: sigma0^2 in column direction (0.0085) |
"SlopeCol":double | Optional parameter, hit variance: diffusion in column direction (0.00011) [0.000397] |
"OffsetDrift":double | Optional parameter, hit variance: sigma0^2 in drift direction (0.25) |
"SlopeDrift":double | Optional parameter, hit variance: diffusion in drift direction (0.0) |
"OffsetRow":double | Optional parameter, hit variance: sigma0^2 in row direction (0.0) |
"SlopeRow":double | Optional parameter, hit variance: diffusion in row direction (0.0) |
"ChargeConversionFactor":double | Optional parameter, charge conversion factor from ADC values to primary electrons (1.) |
"SegmentMatchingChi2Cut":double | Optional parameter, Chi2/Ndf cut for segment matching (Ndf=4 for B off, 5 for B on) (30.) |
"ReferencePointAtPca":bool | Optional parameter, use PCA as reference point, else 1. hit (false) |
"SkipRoadSearch":bool | Optional parameter, skip road search, hit finding only (false) |
"SkipMultiplePulses":bool | Optional parameter, skip multiple pulse candidates as road search seeds (true) |
createProcessor.py
Definition at line 99 of file RowBasedPadPulseRoadSearchProcessor.h.
marlintpc::RowBasedPadPulseRoadSearchProcessor::RowBasedPadPulseRoadSearchProcessor | ( | ) |
Construct processor.
Definition at line 55 of file RowBasedPadPulseRoadSearchProcessor.cc.
References _bfieldScaleFactor, _chargeConversionFactor, _driftVelocity, _dtMax, _inputTrackerPulsesColName, _maxColDiff, _maxColDiffHit, _maxColGap, _maxDxyRoad, _maxDzRoad, _maxOverlapFraction, _maxSeedsInRow, _maxVarXY, _minNumPulses, _minRowDiff, _minRows, _offsetCol, _offsetDrift, _offsetRow, _outputTrackColName, _outputTrackerHitsColName, _refAtPCA, _rowDensityCut, _segCut, _skipMultiplePulses, _skipRoadSearch, _slopeCol, _slopeDrift, and _slopeRow.
Referenced by newProcessor().
|
private |
Check if modules are neighbours.
[in] | mod1 | first module |
[in] | mod2 | second module |
[in] | cartesian | flag for cartesian coordinates |
Use getModuleExtent to define module center and extension, cut on distance normalized to (projected) extension (squared)
Definition at line 882 of file RowBasedPadPulseRoadSearchProcessor.cc.
Referenced by init().
|
virtual |
Definition at line 866 of file RowBasedPadPulseRoadSearchProcessor.cc.
Referenced by newProcessor().
|
virtual |
Definition at line 870 of file RowBasedPadPulseRoadSearchProcessor.cc.
Referenced by newProcessor().
|
virtual |
Initialize processor.
Definition at line 128 of file RowBasedPadPulseRoadSearchProcessor.cc.
References _modNeighbours, _moduleRowOffset, _Xcenter, _Ycenter, and areNeighbourModules().
Referenced by newProcessor().
|
inlinevirtual |
Definition at line 102 of file RowBasedPadPulseRoadSearchProcessor.h.
References check(), end(), init(), processEvent(), processRunHeader(), and RowBasedPadPulseRoadSearchProcessor().
|
virtual |
Process event.
Definition at line 192 of file RowBasedPadPulseRoadSearchProcessor.cc.
References _bfieldScaleFactor, _Bzc, _chargeConversionFactor, _driftVelocity, _dtMax, _inputTrackerPulsesColName, _maxColDiff, _maxColDiffHit, _maxColGap, _maxDxyRoad, _maxDzRoad, _maxOverlapFraction, _maxSeedsInRow, _maxVarXY, _minNumPulses, _minRowDiff, _minRows, _modNeighbours, _moduleRowOffset, _offsetCol, _offsetDrift, _offsetRow, _outputTrackColName, _outputTrackerHitsColName, _refAtPCA, _rowDensityCut, _segCut, _skipMultiplePulses, _skipRoadSearch, _slopeCol, _slopeDrift, _slopeRow, _Xcenter, _Ycenter, marlintpc::simpleEquiClasses::addIndex(), marlintpc::simpleEquiClasses::addMatch(), marlintpc::pb_Pulse::getCharge(), marlintpc::rb_Segment::getChi2(), marlintpc::simpleEquiClasses::getClasses(), marlintpc::pb_Pulse::getColumn(), marlintpc::simpleHelix::getExpectedPlanePos(), marlintpc::rb_Segment::getNdf(), marlintpc::pb_Seed::getPar(), marlintpc::pb_Seed::getRefPoint(), marlintpc::pb_Pulse::getSeedId(), marlintpc::pb_Seed::getSeedId(), marlintpc::pb_Pulse::getTime(), marlintpc::pb_Seed::getValid(), marlintpc::pb_Pulse::inHit(), marlintpc::pulseflag::isAnomalousShape(), marlintpc::pulseflag::isMultiplePulseCandidate(), marlintpc::pulseflag::isOverflow(), marlintpc::pulseflag::isPlateauCutoff(), marlintpc::pulseflag::isSplit(), marlintpc::pb_Seed::markPulses(), marlintpc::pb_Seed::roadSearch(), marlintpc::hitflag::setAnomalousPulseShape(), marlintpc::hitflag::setAtModuleBorder(), marlintpc::hitflag::setDeadChannel(), marlintpc::hitflag::setMultipleHitCandidate(), marlintpc::hitflag::setOverflowPulse(), marlintpc::hitflag::setPlateauCutoffPulse(), and marlintpc::hitflag::setSplitPulse().
Referenced by newProcessor().
|
virtual |
Definition at line 187 of file RowBasedPadPulseRoadSearchProcessor.cc.
Referenced by newProcessor().
|
protected |
scale factor for magnetic field (default: 1.0)
Definition at line 123 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
private |
magnetic field strength (Bz*c)
Definition at line 152 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent().
|
protected |
charge conversion factor from ADC values to primary electrons
Definition at line 144 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
drift velocity
Definition at line 124 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
maximum drift time difference in neighboorhood/hit
Definition at line 126 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
Name of the TPc pulse input collection.
Definition at line 120 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
maximal column difference to pulse in neighbourhood
Definition at line 125 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
maximum distance of pad to seed in hit (if > 0)
Definition at line 135 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
maximum column gap in hit
Definition at line 137 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
maximum interpolation deviation in XY (road search)
Definition at line 129 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
maximum interpolation deviation in Z (road search)
Definition at line 130 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
maximum overlap fraction of segment seeds
Definition at line 134 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
maximum number of seeding pulses in row
Definition at line 127 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
maximum (average) XY variance (unweighted Chi2/ndf) for segment seed
Definition at line 133 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
minimum number of pulses in hit
Definition at line 136 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
minimum row difference for pulse pair seeding road
Definition at line 128 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
minimum number of rows for segment seed
Definition at line 131 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
private |
neighbour modules
Definition at line 153 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by init(), and processEvent().
|
private |
row offsets per module layer
Definition at line 156 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by init(), and processEvent().
|
protected |
measurement variance: sigma0^2 in column direction
Definition at line 138 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
measurement variance: sigma0^2 in drift direction
Definition at line 140 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
measurement variance: sigma0^2 in row direction (use rowHeight^2/12. if < 0.)
Definition at line 142 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
Name of the track output collection.
Definition at line 122 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
Name of the TPC hit output collection.
Definition at line 121 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
Use Pca as reference point (else 1. hit)
Definition at line 146 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
minimum row 'density' for segment seed
Definition at line 132 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
Chi2/Ndf cut for segment matching (Ndf=4 for B off, 5 for B on)
Definition at line 145 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
skip multiple pulse candidates as road search seeds
Definition at line 148 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
skip road search, hit finding only
Definition at line 147 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
measurement variance: diffusion term in column direction
Definition at line 139 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
measurement variance: diffusion term in drift direction
Definition at line 141 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
protected |
measurement variance: diffusion term in row direction
Definition at line 143 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().
|
private |
TPC center, X coordinate.
Definition at line 154 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by init(), and processEvent().
|
private |
TPC center, Y coordinate.
Definition at line 155 of file RowBasedPadPulseRoadSearchProcessor.h.
Referenced by init(), and processEvent().