MyMarlinTPC  170316
Public Member Functions | Protected Attributes | Private Member Functions | Private Attributes | List of all members
marlintpc::RowBasedPadPulseRoadSearchProcessor Class Reference

Simultaneous hit and track finding (pad row based). More...

#include <RowBasedPadPulseRoadSearchProcessor.h>

Inheritance diagram for marlintpc::RowBasedPadPulseRoadSearchProcessor:

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...
 

Detailed Description

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 '[]'.

Parameters
"InputTrackerPulses":stringThe name of the input collection of TPC (pad) pulses ("TPCPulses")
"OutputTrackerHits":stringThe name of the output collection with the found hits ("TPCHits")
"OutputTracks":stringThe name of the output collection with the found tracks ("TPCTracks")
"BFieldScaleFactor":doubleOptional parameter, scales magnetic field (map), use 1.0 (default) for field ON or 0.0 for field OFF (1.) [0.]
"DriftVelocity":doubleOptional parameter, drift velocity (74.)
"MaxColDiffNeighbourhood":intOptional parameter, maximal column difference to pulse in neighbourhood (1) [2]
"MaxTimeDiff":doubleOptional parameter, maximum drift time difference in neighboorhood/hit (80.)
"MaxSeedsInRow":intOptional parameter, maximum number of seeding pulses in row (10)
"MinRowDiff":intOptional parameter, minimum row difference for pulse pair seeding road (5)
"MaxDxyRoad":doubleOptional parameter, maximum interpolation deviation in XY (road search) (2.) [4.]
"MaxDzRoad":doubleOptional parameter, maximum interpolation deviation in Z (road search) (3.)
"MinRowNum":intOptional parameter, minimum number of rows for segment seed (6)
"MinRowDensity":doubleOptional parameter, minimum row 'density' for segment seed (0.8)
"MaxVarXY":doubleOptional parameter, maximum (average) XY variance (unweighted Chi2/ndf) for segment seed (0.25) [2.]
"MaxOverlapFraction":doubleOptional parameter, maximum overlap fraction of segment seeds (0.5)
"MaxColDiffHit":intOptional parameter, maximum distance of pad to seed in hit (if > 0) (0)
"MinNumPulses":doubleOptional parameter, minimum number of pulses in hit (2)
"MaxColGap":doubleOptional parameter, maximum column gap in hit (1)
"OffsetCol":doubleOptional parameter, hit variance: sigma0^2 in column direction (0.0085)
"SlopeCol":doubleOptional parameter, hit variance: diffusion in column direction (0.00011) [0.000397]
"OffsetDrift":doubleOptional parameter, hit variance: sigma0^2 in drift direction (0.25)
"SlopeDrift":doubleOptional parameter, hit variance: diffusion in drift direction (0.0)
"OffsetRow":doubleOptional parameter, hit variance: sigma0^2 in row direction (0.0)
"SlopeRow":doubleOptional parameter, hit variance: diffusion in row direction (0.0)
"ChargeConversionFactor":doubleOptional parameter, charge conversion factor from ADC values to primary electrons (1.)
"SegmentMatchingChi2Cut":doubleOptional parameter, Chi2/Ndf cut for segment matching (Ndf=4 for B off, 5 for B on) (30.)
"ReferencePointAtPca":boolOptional parameter, use PCA as reference point, else 1. hit (false)
"SkipRoadSearch":boolOptional parameter, skip road search, hit finding only (false)
"SkipMultiplePulses":boolOptional parameter, skip multiple pulse candidates as road search seeds (true)
Author
C. Kleinwort (170310)
Credits:
The processor skeleton was generated by the script createProcessor.py
Todo:
for track finding
  • definition of module neighbourhood from POLAR module extend ?
  • better resolve equivalent (indirectly matching) overlapping segments ?
Todo:
for hit finding
  • implement pulse (charge, time and Z) calibration (from database)
  • implement dead and noisy pads (from database)
  • implement more sophisticated resolution of ambiguities and hit splitting

Definition at line 99 of file RowBasedPadPulseRoadSearchProcessor.h.

Constructor & Destructor Documentation

◆ RowBasedPadPulseRoadSearchProcessor()

marlintpc::RowBasedPadPulseRoadSearchProcessor::RowBasedPadPulseRoadSearchProcessor ( )

Member Function Documentation

◆ areNeighbourModules()

bool marlintpc::RowBasedPadPulseRoadSearchProcessor::areNeighbourModules ( gear::TPCModule *  mod1,
gear::TPCModule *  mod2,
bool  cartesian 
)
private

Check if modules are neighbours.

Parameters
[in]mod1first module
[in]mod2second module
[in]cartesianflag for cartesian coordinates
Returns
flag
Todo:
implement non cartesian (=polar) case

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().

◆ check()

void marlintpc::RowBasedPadPulseRoadSearchProcessor::check ( EVENT::LCEvent *  evt)
virtual

Definition at line 866 of file RowBasedPadPulseRoadSearchProcessor.cc.

Referenced by newProcessor().

◆ end()

void marlintpc::RowBasedPadPulseRoadSearchProcessor::end ( )
virtual

Definition at line 870 of file RowBasedPadPulseRoadSearchProcessor.cc.

Referenced by newProcessor().

◆ init()

void marlintpc::RowBasedPadPulseRoadSearchProcessor::init ( )
virtual

Initialize processor.

Definition at line 128 of file RowBasedPadPulseRoadSearchProcessor.cc.

References _modNeighbours, _moduleRowOffset, _Xcenter, _Ycenter, and areNeighbourModules().

Referenced by newProcessor().

◆ newProcessor()

virtual Processor* marlintpc::RowBasedPadPulseRoadSearchProcessor::newProcessor ( )
inlinevirtual

◆ processEvent()

void marlintpc::RowBasedPadPulseRoadSearchProcessor::processEvent ( EVENT::LCEvent *  evt)
virtual

Process event.

  • Prepare input: map of modules with maps of rows with vectors of (pad) pulses.
  • For each module independently build hits and track segments from pad roads:
  1. Define (column) neighbourhoods (for local charge density).
  2. Road search starting with pairs of seeding pulses.
  3. Resolve (large segment) seed overlaps.
  4. Resolve (leading) pulse ambiguities.
  5. Road search for trailing pulses.
  6. Resolve (trailing) pulse ambiguities.
  7. Make hits (in rows seeded by leading pulses)
  8. Make segments (in modules from hits).
  • Link segments between modules.
  • Fill track output collection.

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().

◆ processRunHeader()

void marlintpc::RowBasedPadPulseRoadSearchProcessor::processRunHeader ( EVENT::LCRunHeader *  run)
virtual

Definition at line 187 of file RowBasedPadPulseRoadSearchProcessor.cc.

Referenced by newProcessor().

Member Data Documentation

◆ _bfieldScaleFactor

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_bfieldScaleFactor
protected

scale factor for magnetic field (default: 1.0)

Definition at line 123 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _Bzc

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_Bzc
private

magnetic field strength (Bz*c)

Definition at line 152 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent().

◆ _chargeConversionFactor

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_chargeConversionFactor
protected

charge conversion factor from ADC values to primary electrons

Definition at line 144 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _driftVelocity

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_driftVelocity
protected

drift velocity

Definition at line 124 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _dtMax

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_dtMax
protected

maximum drift time difference in neighboorhood/hit

Definition at line 126 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _inputTrackerPulsesColName

std::string marlintpc::RowBasedPadPulseRoadSearchProcessor::_inputTrackerPulsesColName
protected

Name of the TPc pulse input collection.

Definition at line 120 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _maxColDiff

int marlintpc::RowBasedPadPulseRoadSearchProcessor::_maxColDiff
protected

maximal column difference to pulse in neighbourhood

Definition at line 125 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _maxColDiffHit

int marlintpc::RowBasedPadPulseRoadSearchProcessor::_maxColDiffHit
protected

maximum distance of pad to seed in hit (if > 0)

Definition at line 135 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _maxColGap

int marlintpc::RowBasedPadPulseRoadSearchProcessor::_maxColGap
protected

maximum column gap in hit

Definition at line 137 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _maxDxyRoad

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_maxDxyRoad
protected

maximum interpolation deviation in XY (road search)

Definition at line 129 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _maxDzRoad

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_maxDzRoad
protected

maximum interpolation deviation in Z (road search)

Definition at line 130 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _maxOverlapFraction

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_maxOverlapFraction
protected

maximum overlap fraction of segment seeds

Definition at line 134 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _maxSeedsInRow

int marlintpc::RowBasedPadPulseRoadSearchProcessor::_maxSeedsInRow
protected

maximum number of seeding pulses in row

Definition at line 127 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _maxVarXY

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_maxVarXY
protected

maximum (average) XY variance (unweighted Chi2/ndf) for segment seed

Definition at line 133 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _minNumPulses

int marlintpc::RowBasedPadPulseRoadSearchProcessor::_minNumPulses
protected

minimum number of pulses in hit

Definition at line 136 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _minRowDiff

int marlintpc::RowBasedPadPulseRoadSearchProcessor::_minRowDiff
protected

minimum row difference for pulse pair seeding road

Definition at line 128 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _minRows

int marlintpc::RowBasedPadPulseRoadSearchProcessor::_minRows
protected

minimum number of rows for segment seed

Definition at line 131 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _modNeighbours

std::map<int, std::vector<int> > marlintpc::RowBasedPadPulseRoadSearchProcessor::_modNeighbours
private

neighbour modules

Definition at line 153 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by init(), and processEvent().

◆ _moduleRowOffset

std::map<int, int> marlintpc::RowBasedPadPulseRoadSearchProcessor::_moduleRowOffset
private

row offsets per module layer

Definition at line 156 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by init(), and processEvent().

◆ _offsetCol

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_offsetCol
protected

measurement variance: sigma0^2 in column direction

Definition at line 138 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _offsetDrift

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_offsetDrift
protected

measurement variance: sigma0^2 in drift direction

Definition at line 140 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _offsetRow

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_offsetRow
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().

◆ _outputTrackColName

std::string marlintpc::RowBasedPadPulseRoadSearchProcessor::_outputTrackColName
protected

Name of the track output collection.

Definition at line 122 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _outputTrackerHitsColName

std::string marlintpc::RowBasedPadPulseRoadSearchProcessor::_outputTrackerHitsColName
protected

Name of the TPC hit output collection.

Definition at line 121 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _refAtPCA

bool marlintpc::RowBasedPadPulseRoadSearchProcessor::_refAtPCA
protected

Use Pca as reference point (else 1. hit)

Definition at line 146 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _rowDensityCut

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_rowDensityCut
protected

minimum row 'density' for segment seed

Definition at line 132 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _segCut

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_segCut
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().

◆ _skipMultiplePulses

bool marlintpc::RowBasedPadPulseRoadSearchProcessor::_skipMultiplePulses
protected

skip multiple pulse candidates as road search seeds

Definition at line 148 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _skipRoadSearch

bool marlintpc::RowBasedPadPulseRoadSearchProcessor::_skipRoadSearch
protected

skip road search, hit finding only

Definition at line 147 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _slopeCol

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_slopeCol
protected

measurement variance: diffusion term in column direction

Definition at line 139 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _slopeDrift

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_slopeDrift
protected

measurement variance: diffusion term in drift direction

Definition at line 141 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _slopeRow

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_slopeRow
protected

measurement variance: diffusion term in row direction

Definition at line 143 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by processEvent(), and RowBasedPadPulseRoadSearchProcessor().

◆ _Xcenter

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_Xcenter
private

TPC center, X coordinate.

Definition at line 154 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by init(), and processEvent().

◆ _Ycenter

double marlintpc::RowBasedPadPulseRoadSearchProcessor::_Ycenter
private

TPC center, Y coordinate.

Definition at line 155 of file RowBasedPadPulseRoadSearchProcessor.h.

Referenced by init(), and processEvent().


The documentation for this class was generated from the following files: