CMS 3D CMS Logo

List of all members | Classes | Public Types | Public Member Functions | Private Member Functions | Private Attributes
CSCSegAlgoRU Class Reference

#include <CSCSegAlgoRU.h>

Inheritance diagram for CSCSegAlgoRU:
CSCSegmentAlgorithm

Classes

struct  AlgoState
 

Public Types

typedef std::vector< bool > BoolContainer
 
typedef std::vector< const CSCRecHit2D * > ChamberHitContainer
 
typedef std::vector< const CSCRecHit2D * >::const_iterator ChamberHitContainerCIt
 
typedef std::vector< int > LayerIndex
 
typedef ROOT::Math::SVector< double, 6 > SVector6
 Typedefs. More...
 

Public Member Functions

std::vector< CSCSegmentbuildSegments (const CSCChamber *aChamber, const ChamberHitContainer &rechits) const
 
 CSCSegAlgoRU (const edm::ParameterSet &ps)
 Constructor. More...
 
std::vector< CSCSegmentrun (const CSCChamber *aChamber, const ChamberHitContainer &rechits) override
 
 ~CSCSegAlgoRU () override
 Destructor. More...
 
- Public Member Functions inherited from CSCSegmentAlgorithm
 CSCSegmentAlgorithm (const edm::ParameterSet &)
 Constructor. More...
 
virtual std::vector< CSCSegmentrun (const CSCChamber *chamber, const std::vector< const CSCRecHit2D * > &rechits)=0
 
virtual ~CSCSegmentAlgorithm ()
 Destructor. More...
 

Private Member Functions

bool addHit (AlgoState &aState, const CSCRecHit2D *hit, int layer) const
 Utility functions. More...
 
bool areHitsCloseInGlobalPhi (const AlgoState &aState, const CSCRecHit2D *h1, const CSCRecHit2D *h2) const
 
bool areHitsCloseInR (const AlgoState &aState, const CSCRecHit2D *h1, const CSCRecHit2D *h2) const
 Utility functions. More...
 
void baseline (AlgoState &aState, int n_seg_min) const
 
void compareProtoSegment (AlgoState &aState, const CSCRecHit2D *h, int layer) const
 
float fit_r_phi (const AlgoState &aState, const SVector6 &points, int layer) const
 
float fitX (const AlgoState &aState, SVector6 points, SVector6 errors, int ir, int ir2, float &chi2_str) const
 
void flagHitsAsUsed (const AlgoState &aState, const ChamberHitContainer &rechitsInChamber, BoolContainer &used) const
 
bool hasHitOnLayer (const AlgoState &aState, int layer) const
 
void increaseProtoSegment (AlgoState &aState, const CSCRecHit2D *h, int layer, int chi2_factor) const
 
bool isHitNearSegment (const AlgoState &aState, const CSCRecHit2D *h) const
 
bool isSegmentGood (const AlgoState &aState, const ChamberHitContainer &rechitsInChamber) const
 
float phiAtZ (const AlgoState &aState, float z) const
 
bool replaceHit (AlgoState &aState, const CSCRecHit2D *h, int layer) const
 
void tryAddingHitsToSegment (AlgoState &aState, const ChamberHitContainer &rechitsInChamber, const BoolContainer &used, const LayerIndex &layerIndex, const ChamberHitContainerCIt i1, const ChamberHitContainerCIt i2) const
 
void updateParameters (AlgoState &aState) const
 

Private Attributes

float chi2_str_
 
float chi2Max
 
float chi2Norm_2D_
 
bool debugInfo
 
bool doCollisions
 
float dPhiIntMax
 
float dPhiMax
 
float dRIntMax
 
float dRMax
 
bool enlarge
 
int minLayersApart
 
const std::string myName
 
double theChi2
 
LocalVector theDirection
 
LocalPoint theOrigin
 
float uz
 
float vz
 
float wideSeg
 

Detailed Description

This is the original algorithm for building endcap muon track segments out of the rechit's in a CSCChamber 'RU' = 'RUssia' = Road Usage

A CSCSegment is a RecSegment4D, and is built from CSCRecHit2D objects, each of which is a RecHit2DLocalPos.

This class is used by the CSCSegmentAlgorithm.
Alternative algorithms can be used for the segment building by writing classes like this, and then selecting which one is actually used via the CSCSegmentBuilder.

developed and implemented by Vladimir Palichik Vladi.nosp@m.mir..nosp@m.Paltc.nosp@m.hik@.nosp@m.cern..nosp@m.ch and Nikolay Voytishin nikol.nosp@m.ay.v.nosp@m.oytis.nosp@m.hin@.nosp@m.cern..nosp@m.ch

Definition at line 35 of file CSCSegAlgoRU.h.

Member Typedef Documentation

◆ BoolContainer

typedef std::vector<bool> CSCSegAlgoRU::BoolContainer

Definition at line 59 of file CSCSegAlgoRU.h.

◆ ChamberHitContainer

typedef std::vector<const CSCRecHit2D*> CSCSegAlgoRU::ChamberHitContainer

Definition at line 51 of file CSCSegAlgoRU.h.

◆ ChamberHitContainerCIt

typedef std::vector<const CSCRecHit2D*>::const_iterator CSCSegAlgoRU::ChamberHitContainerCIt

Definition at line 52 of file CSCSegAlgoRU.h.

◆ LayerIndex

typedef std::vector<int> CSCSegAlgoRU::LayerIndex

Definition at line 50 of file CSCSegAlgoRU.h.

◆ SVector6

typedef ROOT::Math::SVector<double, 6> CSCSegAlgoRU::SVector6

Typedefs.

Definition at line 48 of file CSCSegAlgoRU.h.

Constructor & Destructor Documentation

◆ CSCSegAlgoRU()

CSCSegAlgoRU::CSCSegAlgoRU ( const edm::ParameterSet ps)
explicit

Constructor.

Definition at line 22 of file CSCSegAlgoRU.cc.

22  : CSCSegmentAlgorithm(ps), myName("CSCSegAlgoRU") {
23  doCollisions = ps.getParameter<bool>("doCollisions");
24  enlarge = ps.getParameter<bool>("enlarge");
25  chi2_str_ = ps.getParameter<double>("chi2_str");
26  chi2Norm_2D_ = ps.getParameter<double>("chi2Norm_2D_");
27  dRMax = ps.getParameter<double>("dRMax");
28  dPhiMax = ps.getParameter<double>("dPhiMax");
29  dRIntMax = ps.getParameter<double>("dRIntMax");
30  dPhiIntMax = ps.getParameter<double>("dPhiIntMax");
31  chi2Max = ps.getParameter<double>("chi2Max");
32  wideSeg = ps.getParameter<double>("wideSeg");
33  minLayersApart = ps.getParameter<int>("minLayersApart");
34 
35  LogDebug("CSC") << myName << " has algorithm cuts set to: \n"
36  << "--------------------------------------------------------------------\n"
37  << "dRMax = " << dRMax << '\n'
38  << "dPhiMax = " << dPhiMax << '\n'
39  << "dRIntMax = " << dRIntMax << '\n'
40  << "dPhiIntMax = " << dPhiIntMax << '\n'
41  << "chi2Max = " << chi2Max << '\n'
42  << "wideSeg = " << wideSeg << '\n'
43  << "minLayersApart = " << minLayersApart << std::endl;
44 
45  //reset the thresholds for non-collision data
46  if (!doCollisions) {
47  dRMax = 2.0;
48  dPhiMax = 2 * dPhiMax;
49  dRIntMax = 2 * dRIntMax;
50  dPhiIntMax = 2 * dPhiIntMax;
52  chi2_str_ = 100;
53  chi2Max = 2 * chi2Max;
54  }
55 }

References chi2_str_, chi2Max, chi2Norm_2D_, doCollisions, dPhiIntMax, dPhiMax, dRIntMax, dRMax, enlarge, edm::ParameterSet::getParameter(), LogDebug, minLayersApart, myName, and wideSeg.

◆ ~CSCSegAlgoRU()

CSCSegAlgoRU::~CSCSegAlgoRU ( )
inlineoverride

Destructor.

Definition at line 64 of file CSCSegAlgoRU.h.

64 {};

Member Function Documentation

◆ addHit()

bool CSCSegAlgoRU::addHit ( AlgoState aState,
const CSCRecHit2D hit,
int  layer 
) const
private

Utility functions.

Definition at line 608 of file CSCSegAlgoRU.cc.

608  {
609  // Return true if hit was added successfully
610  // (and then parameters are updated).
611  // Return false if there is already a hit on the same layer, or insert failed.
612  ChamberHitContainer::const_iterator it;
613  for (it = aState.proto_segment.begin(); it != aState.proto_segment.end(); it++)
614  if (((*it)->cscDetId().layer() == layer) && (aHit != (*it)))
615  return false;
616  aState.proto_segment.push_back(aHit);
617  // make a fit
618  updateParameters(aState);
619  return true;
620 }

References CSCSegAlgoRU::AlgoState::proto_segment, and updateParameters().

Referenced by buildSegments(), increaseProtoSegment(), and replaceHit().

◆ areHitsCloseInGlobalPhi()

bool CSCSegAlgoRU::areHitsCloseInGlobalPhi ( const AlgoState aState,
const CSCRecHit2D h1,
const CSCRecHit2D h2 
) const
private

Definition at line 444 of file CSCSegAlgoRU.cc.

446  {
447  float strip_width[10] = {0.003878509,
448  0.002958185,
449  0.002327105,
450  0.00152552,
451  0.00465421,
452  0.002327105,
453  0.00465421,
454  0.002327105,
455  0.00465421,
456  0.002327105}; //in rad
457  const CSCLayer* l1 = aState.aChamber->layer(h1->cscDetId().layer());
458  GlobalPoint gp1 = l1->toGlobal(h1->localPosition());
459  const CSCLayer* l2 = aState.aChamber->layer(h2->cscDetId().layer());
460  GlobalPoint gp2 = l2->toGlobal(h2->localPosition());
461  float err_stpos_h1 = h1->errorWithinStrip();
462  float err_stpos_h2 = h2->errorWithinStrip();
463  CSCDetId id = h1->cscDetId();
464  int iStn = id.iChamberType() - 1;
465  float dphi_incr = 0;
466  if (err_stpos_h1 > 0.25 * strip_width[iStn] || err_stpos_h2 > 0.25 * strip_width[iStn])
467  dphi_incr = 0.5 * strip_width[iStn];
468  float dphi12 = deltaPhi(gp1.barePhi(), gp2.barePhi());
469  return (fabs(dphi12) < (aState.dPhiMax * aState.strip_iadd + dphi_incr)) ? true : false; // +v
470 }

References CSCSegAlgoRU::AlgoState::aChamber, PV3DBase< T, PVType, FrameType >::barePhi(), CSCRecHit2D::cscDetId(), SiPixelRawToDigiRegional_cfi::deltaPhi, CSCSegAlgoRU::AlgoState::dPhiMax, CSCRecHit2D::errorWithinStrip(), CSCDetId::iChamberType(), CSCDetId::layer(), CSCChamber::layer(), CSCRecHit2D::localPosition(), CSCSegAlgoRU::AlgoState::strip_iadd, GeomDet::toGlobal(), and funct::true.

Referenced by buildSegments().

◆ areHitsCloseInR()

bool CSCSegAlgoRU::areHitsCloseInR ( const AlgoState aState,
const CSCRecHit2D h1,
const CSCRecHit2D h2 
) const
private

Utility functions.

Definition at line 397 of file CSCSegAlgoRU.cc.

397  {
398  float maxWG_width[10] = {0, 0, 4.1, 5.69, 2.49, 5.06, 2.49, 5.06, 1.87, 5.06};
399  CSCDetId id = h1->cscDetId();
400  int iStn = id.iChamberType() - 1;
401  //find maxWG_width for ME11 (tilt = 29deg)
402  int wg_num = h2->hitWire();
403  if (iStn == 0 || iStn == 1) {
404  if (wg_num == 1) {
405  maxWG_width[0] = 9.25;
406  maxWG_width[1] = 9.25;
407  }
408  if (wg_num > 1 && wg_num < 48) {
409  maxWG_width[0] = 3.14;
410  maxWG_width[1] = 3.14;
411  }
412  if (wg_num == 48) {
413  maxWG_width[0] = 10.75;
414  maxWG_width[1] = 10.75;
415  }
416  }
417  const CSCLayer* l1 = aState.aChamber->layer(h1->cscDetId().layer());
418  GlobalPoint gp1 = l1->toGlobal(h1->localPosition());
419  const CSCLayer* l2 = aState.aChamber->layer(h2->cscDetId().layer());
420  GlobalPoint gp2 = l2->toGlobal(h2->localPosition());
421  //find z to understand the direction
422  float h1z = gp1.z();
423  float h2z = gp2.z();
424  //switch off the IP check for non collisions case
425  if (!aState.doCollisions) {
426  h1z = 1;
427  h2z = 1;
428  }
429 
430  if (aState.enlarge) {
431  return (gp2.perp() > ((gp1.perp() - aState.dRMax * aState.strip_iadd * maxWG_width[iStn]) * h2z) / h1z &&
432  gp2.perp() < ((gp1.perp() + aState.dRMax * aState.strip_iadd * maxWG_width[iStn]) * h2z) / h1z)
433  ? true
434  : false;
435 
436  } else {
437  return (gp2.perp() > ((gp1.perp() - aState.dRMax * maxWG_width[iStn]) * h2z) / h1z &&
438  gp2.perp() < ((gp1.perp() + aState.dRMax * maxWG_width[iStn]) * h2z) / h1z)
439  ? true
440  : false;
441  }
442 }

References CSCSegAlgoRU::AlgoState::aChamber, CSCRecHit2D::cscDetId(), CSCSegAlgoRU::AlgoState::doCollisions, CSCSegAlgoRU::AlgoState::dRMax, CSCSegAlgoRU::AlgoState::enlarge, CSCRecHit2D::hitWire(), CSCDetId::iChamberType(), CSCDetId::layer(), CSCChamber::layer(), CSCRecHit2D::localPosition(), PV3DBase< T, PVType, FrameType >::perp(), CSCSegAlgoRU::AlgoState::strip_iadd, GeomDet::toGlobal(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by buildSegments().

◆ baseline()

void CSCSegAlgoRU::baseline ( AlgoState aState,
int  n_seg_min 
) const
private

(nhits-2)

(nhits-3)

Definition at line 649 of file CSCSegAlgoRU.cc.

649  {
650  int nhits = aState.proto_segment.size();
651  //initialise vectors for strip position and error within strip
652  SVector6 sp;
653  SVector6 se;
654  unsigned int init_size = aState.proto_segment.size();
656  buffer.clear();
657  buffer.reserve(init_size);
658  while (buffer.size() < init_size) {
659  ChamberHitContainer::iterator min;
660  int min_layer = 10;
661  for (ChamberHitContainer::iterator k = aState.proto_segment.begin(); k != aState.proto_segment.end(); k++) {
662  const CSCRecHit2D* iRHk = *k;
663  CSCDetId idRHk = iRHk->cscDetId();
664  int kLayer = idRHk.layer();
665  if (kLayer < min_layer) {
666  min_layer = kLayer;
667  min = k;
668  }
669  }
670  buffer.push_back(*min);
671  aState.proto_segment.erase(min);
672  } //while
673 
674  aState.proto_segment.clear();
675  for (ChamberHitContainer::const_iterator cand = buffer.begin(); cand != buffer.end(); cand++) {
676  aState.proto_segment.push_back(*cand);
677  }
678 
679  for (ChamberHitContainer::const_iterator iRH = aState.proto_segment.begin(); iRH != aState.proto_segment.end();
680  iRH++) {
681  const CSCRecHit2D* iRHp = *iRH;
682  CSCDetId idRH = iRHp->cscDetId();
683  int kRing = idRH.ring();
684  int kStation = idRH.station();
685  int kLayer = idRH.layer();
686  // Find the strip containing this hit
687  int centerid = iRHp->nStrips() / 2;
688  int centerStrip = iRHp->channels(centerid);
689  float stpos = (*iRHp).positionWithinStrip();
690  se(kLayer - 1) = (*iRHp).errorWithinStrip();
691  // Take into account half-strip staggering of layers (ME1/1 has no staggering)
692  if (kStation == 1 && (kRing == 1 || kRing == 4))
693  sp(kLayer - 1) = stpos + centerStrip;
694  else {
695  if (kLayer == 1 || kLayer == 3 || kLayer == 5)
696  sp(kLayer - 1) = stpos + centerStrip;
697  if (kLayer == 2 || kLayer == 4 || kLayer == 6)
698  sp(kLayer - 1) = stpos - 0.5 + centerStrip;
699  }
700  }
701  float chi2_str;
702  fitX(aState, sp, se, -1, -1, chi2_str);
703 
704  //-----------------------------------------------------
705  // Optimal point rejection method
706  //-----------------------------------------------------
707  float minSum = 1000;
708  int i1b = 0;
709  int i2b = 0;
710  int iworst = -1;
711  int bad_layer = -1;
712  ChamberHitContainer::const_iterator rh_to_be_deleted_1;
713  ChamberHitContainer::const_iterator rh_to_be_deleted_2;
714  if ((chi2_str) > aState.chi2_str_ * aState.chi2D_iadd) {
715  for (ChamberHitContainer::const_iterator i1 = aState.proto_segment.begin(); i1 != aState.proto_segment.end();
716  ++i1) {
717  ++i1b;
718  const CSCRecHit2D* i1_1 = *i1;
719  CSCDetId idRH1 = i1_1->cscDetId();
720  int z1 = idRH1.layer();
721  i2b = i1b;
722  for (ChamberHitContainer::const_iterator i2 = i1 + 1; i2 != aState.proto_segment.end(); ++i2) {
723  ++i2b;
724  const CSCRecHit2D* i2_1 = *i2;
725  CSCDetId idRH2 = i2_1->cscDetId();
726  int z2 = idRH2.layer();
727  int irej = 0;
728  for (ChamberHitContainer::const_iterator ir = aState.proto_segment.begin(); ir != aState.proto_segment.end();
729  ++ir) {
730  ++irej;
731  if (ir == i1 || ir == i2)
732  continue;
733  float dsum = 0;
734  int hit_nr = 0;
735  const CSCRecHit2D* ir_1 = *ir;
736  CSCDetId idRH = ir_1->cscDetId();
737  int worst_layer = idRH.layer();
738  for (ChamberHitContainer::const_iterator i = aState.proto_segment.begin(); i != aState.proto_segment.end();
739  ++i) {
740  ++hit_nr;
741  const CSCRecHit2D* i_1 = *i;
742  if (i == i1 || i == i2 || i == ir)
743  continue;
744  float slope = (sp(z2 - 1) - sp(z1 - 1)) / (z2 - z1);
745  float intersept = sp(z1 - 1) - slope * z1;
746  CSCDetId idRH = i_1->cscDetId();
747  int z = idRH.layer();
748  float di = fabs(sp(z - 1) - intersept - slope * z);
749  dsum = dsum + di;
750  } //i
751  if (dsum < minSum) {
752  minSum = dsum;
753  bad_layer = worst_layer;
754  iworst = irej;
755  rh_to_be_deleted_1 = ir;
756  }
757  } //ir
758  } //i2
759  } //i1
760  fitX(aState, sp, se, bad_layer, -1, chi2_str);
761  } //if chi2prob<1.0e-4
762 
763  //find worst from n-1 hits
764  int iworst2 = -1;
765  int bad_layer2 = -1;
766  if (iworst > -1 && (nhits - 1) > n_seg_min && (chi2_str) > aState.chi2_str_ * aState.chi2D_iadd) {
767  iworst = -1;
768  float minSum = 1000;
769  int i1b = 0;
770  int i2b = 0;
771  for (ChamberHitContainer::const_iterator i1 = aState.proto_segment.begin(); i1 != aState.proto_segment.end();
772  ++i1) {
773  ++i1b;
774  const CSCRecHit2D* i1_1 = *i1;
775  CSCDetId idRH1 = i1_1->cscDetId();
776  int z1 = idRH1.layer();
777  i2b = i1b;
778  for (ChamberHitContainer::const_iterator i2 = i1 + 1; i2 != aState.proto_segment.end(); ++i2) {
779  ++i2b;
780  const CSCRecHit2D* i2_1 = *i2;
781  CSCDetId idRH2 = i2_1->cscDetId();
782  int z2 = idRH2.layer();
783  int irej = 0;
784  for (ChamberHitContainer::const_iterator ir = aState.proto_segment.begin(); ir != aState.proto_segment.end();
785  ++ir) {
786  ++irej;
787  int irej2 = 0;
788  if (ir == i1 || ir == i2)
789  continue;
790  const CSCRecHit2D* ir_1 = *ir;
791  CSCDetId idRH = ir_1->cscDetId();
792  int worst_layer = idRH.layer();
793  for (ChamberHitContainer::const_iterator ir2 = aState.proto_segment.begin();
794  ir2 != aState.proto_segment.end();
795  ++ir2) {
796  ++irej2;
797  if (ir2 == i1 || ir2 == i2 || ir2 == ir)
798  continue;
799  float dsum = 0;
800  int hit_nr = 0;
801  const CSCRecHit2D* ir2_1 = *ir2;
802  CSCDetId idRH = ir2_1->cscDetId();
803  int worst_layer2 = idRH.layer();
804  for (ChamberHitContainer::const_iterator i = aState.proto_segment.begin(); i != aState.proto_segment.end();
805  ++i) {
806  ++hit_nr;
807  const CSCRecHit2D* i_1 = *i;
808  if (i == i1 || i == i2 || i == ir || i == ir2)
809  continue;
810  float slope = (sp(z2 - 1) - sp(z1 - 1)) / (z2 - z1);
811  float intersept = sp(z1 - 1) - slope * z1;
812  CSCDetId idRH = i_1->cscDetId();
813  int z = idRH.layer();
814  float di = fabs(sp(z - 1) - intersept - slope * z);
815  dsum = dsum + di;
816  } //i
817  if (dsum < minSum) {
818  minSum = dsum;
819  iworst2 = irej2;
820  iworst = irej;
821  bad_layer = worst_layer;
822  bad_layer2 = worst_layer2;
823  rh_to_be_deleted_1 = ir;
824  rh_to_be_deleted_2 = ir2;
825  }
826  } //ir2
827  } //ir
828  } //i2
829  } //i1
830  fitX(aState, sp, se, bad_layer, bad_layer2, chi2_str);
831  } //if prob(n-1)<e-4
832 
833  //----------------------------------
834  //erase bad_hits
835  //----------------------------------
836  if (iworst2 - 1 >= 0 && iworst2 <= int(aState.proto_segment.size())) {
837  aState.proto_segment.erase(rh_to_be_deleted_2);
838  }
839  if (iworst - 1 >= 0 && iworst <= int(aState.proto_segment.size())) {
840  aState.proto_segment.erase(rh_to_be_deleted_1);
841  }
842 }

References edmScanValgrind::buffer, CSCRecHit2D::channels(), CSCSegmentAlgorithmRU_cfi::chi2_str, CSCSegAlgoRU::AlgoState::chi2_str_, CSCSegAlgoRU::AlgoState::chi2D_iadd, CSCRecHit2D::cscDetId(), fitX(), mps_fire::i, testProducerWithPsetDescEmpty_cfi::i1, testProducerWithPsetDescEmpty_cfi::i2, dqmdumpme::k, kLayer(), CSCDetId::layer(), min(), nhits, CSCRecHit2D::nStrips(), CSCSegAlgoRU::AlgoState::proto_segment, CSCDetId::ring(), slope, CSCDetId::station(), z, and testProducerWithPsetDescEmpty_cfi::z2.

Referenced by buildSegments().

◆ buildSegments()

std::vector< CSCSegment > CSCSegAlgoRU::buildSegments ( const CSCChamber aChamber,
const ChamberHitContainer rechits 
) const

Build track segments in this chamber (this is where the actual segment-building algorithm hides.)

Definition at line 57 of file CSCSegAlgoRU.cc.

58  {
59  ChamberHitContainer rechits = urechits;
60  LayerIndex layerIndex(rechits.size());
61  int recHits_per_layer[6] = {0, 0, 0, 0, 0, 0};
62  //skip events with high multiplicity of hits
63  if (rechits.size() > 150) {
64  return std::vector<CSCSegment>();
65  }
66  int iadd = 0;
67  for (unsigned int i = 0; i < rechits.size(); i++) {
68  recHits_per_layer[rechits[i]->cscDetId().layer() - 1]++; //count rh per chamber
69  layerIndex[i] = rechits[i]->cscDetId().layer();
70  }
71  double z1 = aChamber->layer(1)->position().z();
72  double z6 = aChamber->layer(6)->position().z();
73  if (std::abs(z1) > std::abs(z6)) {
74  reverse(layerIndex.begin(), layerIndex.end());
75  reverse(rechits.begin(), rechits.end());
76  }
77  if (rechits.size() < 2) {
78  return std::vector<CSCSegment>();
79  }
80  // We have at least 2 hits. We intend to try all possible pairs of hits to start
81  // segment building. 'All possible' means each hit lies on different layers in the chamber.
82  // after all same size segs are build we get rid of the overcrossed segments using the chi2 criteria
83  // the hits from the segs that are left are marked as used and are not added to segs in future iterations
84  // the hits from 3p segs are marked as used separately in order to try to assamble them in longer segments
85  // in case there is a second pass
86  // Choose first hit (as close to IP as possible) h1 and second hit
87  // (as far from IP as possible) h2 To do this we iterate over hits
88  // in the chamber by layer - pick two layers. Then we
89  // iterate over hits within each of these layers and pick h1 and h2
90  // these. If they are 'close enough' we build an empty
91  // segment. Then try adding hits to this segment.
92  // Initialize flags that a given hit has been allocated to a segment
93  BoolContainer used(rechits.size(), false);
94  BoolContainer used3p(rechits.size(), false);
95  // This is going to point to fits to hits, and its content will be used to create a CSCSegment
96  AlgoState aState;
97  aState.aChamber = aChamber;
98  aState.doCollisions = doCollisions;
99  aState.enlarge = enlarge;
100  aState.dRMax = dRMax;
101  aState.dPhiMax = dPhiMax;
102  aState.dRIntMax = dRIntMax;
103  aState.dPhiIntMax = dPhiIntMax;
104  aState.chi2Norm_2D_ = chi2Norm_2D_;
105  aState.chi2_str_ = chi2_str_;
106  aState.chi2Max = chi2Max;
107 
108  int scale_factor = 1;
109  if (aState.enlarge)
110  scale_factor = 2;
111 
112  // Define buffer for segments we build
113  std::vector<CSCSegment> segments;
115  ChamberHitContainerCIt ie = rechits.end();
116  // Possibly allow 3 passes, second widening scale factor for cuts, third for segments from displaced vertices
117  aState.windowScale = 1.; // scale factor for cuts
118  bool search_disp = false;
119  aState.strip_iadd = 1;
120  aState.chi2D_iadd = 1;
121  int npass = (wideSeg > 1.) ? 3 : 2;
122  for (int ipass = 0; ipass < npass; ++ipass) {
123  if (aState.windowScale > 1.) {
124  iadd = 1;
125  aState.strip_iadd = 2 * scale_factor;
126  aState.chi2D_iadd = 2 * scale_factor;
127  if (aState.enlarge) {
128  aState.chi2Max = 2 * chi2Max;
129  if (rechits.size() <= 12)
130  iadd = 0; //allow 3 hit segments for low hit multiplicity chambers
131  }
132  }
133 
134  int used_rh = 0;
135  for (ChamberHitContainerCIt i1 = ib; i1 != ie; ++i1) {
136  if (used[i1 - ib])
137  used_rh++;
138  }
139 
140  //change the tresholds if it's time to look for displaced mu segments
141  if (aState.doCollisions && search_disp &&
142  int(rechits.size() - used_rh) >
143  2) { //check if there are enough recHits left to build a segment from displaced vertices
144  aState.doCollisions = false;
145  aState.windowScale = 1.; // scale factor for cuts
146  aState.dRMax = scale_factor * 2.0;
147  aState.dPhiMax = scale_factor * 2 * aState.dPhiMax;
148  aState.dRIntMax = scale_factor * 2 * aState.dRIntMax;
149  aState.dPhiIntMax = scale_factor * 2 * aState.dPhiIntMax;
150  aState.chi2Norm_2D_ = scale_factor * 5 * aState.chi2Norm_2D_;
151  aState.chi2_str_ = scale_factor * 100;
152  aState.chi2Max = scale_factor * 2 * aState.chi2Max;
153  } else {
154  search_disp = false; //make sure the flag is off
155  }
156 
157  for (unsigned int n_seg_min = 6u; n_seg_min > 2u + iadd; --n_seg_min) {
158  BoolContainer common_used(rechits.size(), false);
159  std::array<BoolContainer, 120> common_used_it = {};
160  for (unsigned int i = 0; i < common_used_it.size(); i++) {
161  common_used_it[i] = common_used;
162  }
163  ChamberHitContainer best_proto_segment[120];
164  float min_chi[120] = {9999};
165  int common_it = 0;
166  bool first_proto_segment = true;
167  for (ChamberHitContainerCIt i1 = ib; i1 != ie; ++i1) {
168  bool segok = false;
169  //skip if rh is used and the layer tat has big rh multiplicity(>25RHs)
170  if (used[i1 - ib] || recHits_per_layer[int(layerIndex[i1 - ib]) - 1] > 25 ||
171  (n_seg_min == 3 && used3p[i1 - ib]))
172  continue;
173  int layer1 = layerIndex[i1 - ib];
174  const CSCRecHit2D* h1 = *i1;
175  for (ChamberHitContainerCIt i2 = ie - 1; i2 != i1; --i2) {
176  if (used[i2 - ib] || recHits_per_layer[int(layerIndex[i2 - ib]) - 1] > 25 ||
177  (n_seg_min == 3 && used3p[i2 - ib]))
178  continue;
179  int layer2 = layerIndex[i2 - ib];
180  if ((abs(layer2 - layer1) + 1) < int(n_seg_min))
181  break; //decrease n_seg_min
182  const CSCRecHit2D* h2 = *i2;
183  if (areHitsCloseInR(aState, h1, h2) && areHitsCloseInGlobalPhi(aState, h1, h2)) {
184  aState.proto_segment.clear();
185  if (!addHit(aState, h1, layer1))
186  continue;
187  if (!addHit(aState, h2, layer2))
188  continue;
189  // Can only add hits if already have a segment
190  if (aState.sfit)
191  tryAddingHitsToSegment(aState, rechits, used, layerIndex, i1, i2);
192  segok = isSegmentGood(aState, rechits);
193  if (segok) {
194  if (aState.proto_segment.size() > n_seg_min) {
195  baseline(aState, n_seg_min);
196  updateParameters(aState);
197  }
198  if (aState.sfit->chi2() > aState.chi2Norm_2D_ * aState.chi2D_iadd ||
199  aState.proto_segment.size() < n_seg_min)
200  aState.proto_segment.clear();
201  if (!aState.proto_segment.empty()) {
202  updateParameters(aState);
203  //add same-size overcrossed protosegments to the collection
204  if (first_proto_segment) {
205  flagHitsAsUsed(aState, rechits, common_used_it[0]);
206  min_chi[0] = aState.sfit->chi2();
207  best_proto_segment[0] = aState.proto_segment;
208  first_proto_segment = false;
209  } else { //for the rest of found proto_segments
210  common_it++;
211  flagHitsAsUsed(aState, rechits, common_used_it[common_it]);
212  min_chi[common_it] = aState.sfit->chi2();
213  best_proto_segment[common_it] = aState.proto_segment;
214  ChamberHitContainerCIt hi, iu, ik;
215  int iter = common_it;
216  for (iu = ib; iu != ie; ++iu) {
217  for (hi = aState.proto_segment.begin(); hi != aState.proto_segment.end(); ++hi) {
218  if (*hi == *iu) {
219  int merge_nr = -1;
220  for (int k = 0; k < iter + 1; k++) {
221  if (common_used_it[k][iu - ib] == true) {
222  if (merge_nr != -1) {
223  //merge the k and merge_nr collections of flaged hits into the merge_nr collection and unmark the k collection hits
224  for (ik = ib; ik != ie; ++ik) {
225  if (common_used_it[k][ik - ib] == true) {
226  common_used_it[merge_nr][ik - ib] = true;
227  common_used_it[k][ik - ib] = false;
228  }
229  }
230  //change best_protoseg if min_chi_2 is smaller
231  if (min_chi[k] < min_chi[merge_nr]) {
232  min_chi[merge_nr] = min_chi[k];
233  best_proto_segment[merge_nr] = best_proto_segment[k];
234  best_proto_segment[k].clear();
235  min_chi[k] = 9999;
236  }
237  common_it--;
238  } else {
239  merge_nr = k;
240  }
241  } //if(common_used[k][iu-ib] == true)
242  } //for k
243  } //if
244  } //for proto_seg
245  } //for rec_hits
246  } //else
247  } //proto seg not empty
248  }
249  } // h1 & h2 close
250  if (segok)
251  break;
252  } // i2
253  } // i1
254 
255  //add the reconstructed segments
256  for (int j = 0; j < common_it + 1; j++) {
257  aState.proto_segment = best_proto_segment[j];
258  best_proto_segment[j].clear();
259  //SKIP empty proto-segments
260  if (aState.proto_segment.empty())
261  continue;
262  updateParameters(aState);
263  // Create an actual CSCSegment - retrieve all info from the fit
264  CSCSegment temp(aState.sfit->hits(),
265  aState.sfit->intercept(),
266  aState.sfit->localdir(),
267  aState.sfit->covarianceMatrix(),
268  aState.sfit->chi2());
269  aState.sfit = nullptr;
270  segments.push_back(temp);
271  //if the segment has 3 hits flag them as used in a particular way
272  if (aState.proto_segment.size() == 3) {
273  flagHitsAsUsed(aState, rechits, used3p);
274  } else {
275  flagHitsAsUsed(aState, rechits, used);
276  }
277  aState.proto_segment.clear();
278  }
279  } //for n_seg_min
280 
281  if (search_disp) {
282  //reset params and flags for the next chamber
283  search_disp = false;
284  aState.doCollisions = true;
285  aState.dRMax = 2.0;
286  aState.chi2_str_ = 100;
287  aState.dPhiMax = 0.5 * aState.dPhiMax / scale_factor;
288  aState.dRIntMax = 0.5 * aState.dRIntMax / scale_factor;
289  aState.dPhiIntMax = 0.5 * aState.dPhiIntMax / scale_factor;
290  aState.chi2Norm_2D_ = 0.2 * aState.chi2Norm_2D_ / scale_factor;
291  aState.chi2Max = 0.5 * aState.chi2Max / scale_factor;
292  }
293 
294  std::vector<CSCSegment>::iterator it = segments.begin();
295  bool good_segs = false;
296  while (it != segments.end()) {
297  if ((*it).nRecHits() > 3) {
298  good_segs = true;
299  break;
300  }
301  ++it;
302  }
303  if (good_segs &&
304  aState.doCollisions) { // only change window if not enough good segments were found (bool can be changed to int if a >0 number of good segs is required)
305  search_disp = true;
306  continue; //proceed to search the segs from displaced vertices
307  }
308 
309  // Increase cut windows by factor of wideSeg only for collisions
310  if (!aState.doCollisions && !search_disp)
311  break;
312  aState.windowScale = wideSeg;
313  } // ipass
314 
315  //get rid of enchansed 3p segments
316  std::vector<CSCSegment>::iterator it = segments.begin();
317  while (it != segments.end()) {
318  if ((*it).nRecHits() == 3) {
319  bool found_common = false;
320  const std::vector<CSCRecHit2D>& theseRH = (*it).specificRecHits();
321  for (ChamberHitContainerCIt i1 = ib; i1 != ie; ++i1) {
322  if (used[i1 - ib] && used3p[i1 - ib]) {
323  const CSCRecHit2D* sh1 = *i1;
324  CSCDetId id = sh1->cscDetId();
325  int sh1layer = id.layer();
326  int RH_centerid = sh1->nStrips() / 2;
327  int RH_centerStrip = sh1->channels(RH_centerid);
328  int RH_wg = sh1->hitWire();
329  std::vector<CSCRecHit2D>::const_iterator sh;
330  for (sh = theseRH.begin(); sh != theseRH.end(); ++sh) {
331  CSCDetId idRH = sh->cscDetId();
332  //find segment hit coord
333  int shlayer = idRH.layer();
334  int SegRH_centerid = sh->nStrips() / 2;
335  int SegRH_centerStrip = sh->channels(SegRH_centerid);
336  int SegRH_wg = sh->hitWire();
337  if (sh1layer == shlayer && SegRH_centerStrip == RH_centerStrip && SegRH_wg == RH_wg) {
338  //remove the enchansed 3p segment
339  segments.erase(it, (it + 1));
340  found_common = true;
341  break;
342  }
343  } //theserh
344  }
345  if (found_common)
346  break; //current seg has already been erased
347  } //camber hits
348  if (!found_common)
349  ++it;
350  } //its a 3p seg
351  else {
352  ++it; //go to the next seg
353  }
354  } //while
355  // Give the segments to the CSCChamber
356  return segments;
357 } //build segments

References funct::abs(), CSCSegAlgoRU::AlgoState::aChamber, addHit(), areHitsCloseInGlobalPhi(), areHitsCloseInR(), baseline(), CSCRecHit2D::channels(), chi2_str_, chi2Max, chi2Norm_2D_, CSCRecHit2D::cscDetId(), doCollisions, dPhiIntMax, dPhiMax, dRIntMax, dRMax, enlarge, flagHitsAsUsed(), CSCRecHit2D::hitWire(), mps_fire::i, testProducerWithPsetDescEmpty_cfi::i1, testProducerWithPsetDescEmpty_cfi::i2, cuy::ib, isSegmentGood(), dqmiolumiharvest::j, dqmdumpme::k, CSCDetId::layer(), CSCChamber::layer(), CSCRecHit2D::nStrips(), GeomDet::position(), HI_PhotonSkim_cff::rechits, groupFilesInBlocks::reverse, groupFilesInBlocks::temp, tryAddingHitsToSegment(), updateParameters(), wideSeg, and PV3DBase< T, PVType, FrameType >::z().

Referenced by run().

◆ compareProtoSegment()

void CSCSegAlgoRU::compareProtoSegment ( AlgoState aState,
const CSCRecHit2D h,
int  layer 
) const
private

Definition at line 899 of file CSCSegAlgoRU.cc.

899  {
900  // Copy the input CSCSegFit
901  std::unique_ptr<CSCSegFit> oldfit;
902  oldfit = std::make_unique<CSCSegFit>(aState.aChamber, aState.proto_segment);
903  oldfit->fit();
904  ChamberHitContainer oldproto = aState.proto_segment;
905 
906  // May create a new fit
907  bool ok = replaceHit(aState, h, layer);
908  if ((aState.sfit->chi2() >= oldfit->chi2()) || !ok) {
909  // keep original fit
910  aState.proto_segment = oldproto;
911  aState.sfit = std::move(oldfit); // reset to the original input fit
912  }
913 }

References CSCSegAlgoRU::AlgoState::aChamber, eostools::move(), convertSQLiteXML::ok, CSCSegAlgoRU::AlgoState::proto_segment, replaceHit(), and CSCSegAlgoRU::AlgoState::sfit.

Referenced by tryAddingHitsToSegment().

◆ fit_r_phi()

float CSCSegAlgoRU::fit_r_phi ( const AlgoState aState,
const SVector6 points,
int  layer 
) const
private

Definition at line 629 of file CSCSegAlgoRU.cc.

629  {
630  //find R or Phi on the given layer using the given points for the interpolation
631  float Sx = 0;
632  float Sy = 0;
633  float Sxx = 0;
634  float Sxy = 0;
635  for (int i = 1; i < 7; i++) {
636  if (points(i - 1) == 0.)
637  continue;
638  Sy = Sy + (points(i - 1));
639  Sx = Sx + i;
640  Sxx = Sxx + (i * i);
641  Sxy = Sxy + ((i)*points(i - 1));
642  }
643  float delta = 2 * Sxx - Sx * Sx;
644  float intercept = (Sxx * Sy - Sx * Sxy) / delta;
645  float slope = (2 * Sxy - Sx * Sy) / delta;
646  return (intercept + slope * layer);
647 }

References dumpMFGeometry_cfg::delta, mps_fire::i, HLT_FULL_cff::points, and slope.

Referenced by isHitNearSegment().

◆ fitX()

float CSCSegAlgoRU::fitX ( const AlgoState aState,
SVector6  points,
SVector6  errors,
int  ir,
int  ir2,
float &  chi2_str 
) const
private

Definition at line 844 of file CSCSegAlgoRU.cc.

845  {
846  float S = 0;
847  float Sx = 0;
848  float Sy = 0;
849  float Sxx = 0;
850  float Sxy = 0;
851  float sigma2 = 0;
852  for (int i = 1; i < 7; i++) {
853  if (i == ir || i == ir2 || points(i - 1) == 0.)
854  continue;
855  sigma2 = errors(i - 1) * errors(i - 1);
856  float i1 = i - 3.5;
857  S = S + (1 / sigma2);
858  Sy = Sy + (points(i - 1) / sigma2);
859  Sx = Sx + ((i1) / sigma2);
860  Sxx = Sxx + (i1 * i1) / sigma2;
861  Sxy = Sxy + (((i1)*points(i - 1)) / sigma2);
862  }
863  float delta = S * Sxx - Sx * Sx;
864  float intercept = (Sxx * Sy - Sx * Sxy) / delta;
865  float slope = (S * Sxy - Sx * Sy) / delta;
866  float chi_str = 0;
867  chi2_str = 0;
868  // calculate chi2_str
869  for (int i = 1; i < 7; i++) {
870  if (i == ir || i == ir2 || points(i - 1) == 0.)
871  continue;
872  chi_str = (points(i - 1) - intercept - slope * (i - 3.5)) / (errors(i - 1));
873  chi2_str = chi2_str + chi_str * chi_str;
874  }
875  return (intercept + slope * 0);
876 }

References CSCSegmentAlgorithmRU_cfi::chi2_str, dumpMFGeometry_cfg::delta, debug_messages_cfi::errors, mps_fire::i, testProducerWithPsetDescEmpty_cfi::i1, HLT_FULL_cff::points, and slope.

Referenced by baseline().

◆ flagHitsAsUsed()

void CSCSegAlgoRU::flagHitsAsUsed ( const AlgoState aState,
const ChamberHitContainer rechitsInChamber,
BoolContainer used 
) const
private

Flag hits on segment as used

Definition at line 594 of file CSCSegAlgoRU.cc.

596  {
597  // Flag hits on segment as used
598  ChamberHitContainerCIt ib = rechitsInChamber.begin();
600  for (hi = aState.proto_segment.begin(); hi != aState.proto_segment.end(); ++hi) {
601  for (iu = ib; iu != rechitsInChamber.end(); ++iu) {
602  if (*hi == *iu)
603  used[iu - ib] = true;
604  }
605  }
606 }

References cuy::ib, and CSCSegAlgoRU::AlgoState::proto_segment.

Referenced by buildSegments().

◆ hasHitOnLayer()

bool CSCSegAlgoRU::hasHitOnLayer ( const AlgoState aState,
int  layer 
) const
private

Definition at line 878 of file CSCSegAlgoRU.cc.

878  {
879  // Is there is already a hit on this layer?
881  for (it = aState.proto_segment.begin(); it != aState.proto_segment.end(); it++)
882  if ((*it)->cscDetId().layer() == layer)
883  return true;
884  return false;
885 }

References CSCSegAlgoRU::AlgoState::proto_segment.

Referenced by tryAddingHitsToSegment().

◆ increaseProtoSegment()

void CSCSegAlgoRU::increaseProtoSegment ( AlgoState aState,
const CSCRecHit2D h,
int  layer,
int  chi2_factor 
) const
private

Definition at line 915 of file CSCSegAlgoRU.cc.

915  {
916  // Creates a new fit
917  std::unique_ptr<CSCSegFit> oldfit;
918  ChamberHitContainer oldproto = aState.proto_segment;
919  oldfit = std::make_unique<CSCSegFit>(aState.aChamber, aState.proto_segment);
920  oldfit->fit();
921 
922  bool ok = addHit(aState, h, layer);
923  //@@ TEST ON ndof<=0 IS JUST TO ACCEPT nhits=2 CASE??
924  if (!ok || ((aState.sfit->ndof() > 0) && (aState.sfit->chi2() / aState.sfit->ndof() >= aState.chi2Max))) {
925  // reset to original fit
926  aState.proto_segment = oldproto;
927  aState.sfit = std::move(oldfit);
928  }
929 }

References CSCSegAlgoRU::AlgoState::aChamber, addHit(), CSCSegAlgoRU::AlgoState::chi2Max, eostools::move(), convertSQLiteXML::ok, CSCSegAlgoRU::AlgoState::proto_segment, and CSCSegAlgoRU::AlgoState::sfit.

Referenced by tryAddingHitsToSegment().

◆ isHitNearSegment()

bool CSCSegAlgoRU::isHitNearSegment ( const AlgoState aState,
const CSCRecHit2D h 
) const
private

Definition at line 472 of file CSCSegAlgoRU.cc.

472  {
473  // Is hit near segment?
474  // Requires deltaphi and deltaR within ranges specified in parameter set.
475  // Note that to make intuitive cuts on delta(phi) one must work in
476  // phi range (-pi, +pi] not [0, 2pi)
477  float strip_width[10] = {0.003878509,
478  0.002958185,
479  0.002327105,
480  0.00152552,
481  0.00465421,
482  0.002327105,
483  0.00465421,
484  0.002327105,
485  0.00465421,
486  0.002327105}; //in rad
487  const CSCLayer* l1 = aState.aChamber->layer((*(aState.proto_segment.begin()))->cscDetId().layer());
488  GlobalPoint gp1 = l1->toGlobal((*(aState.proto_segment.begin()))->localPosition());
489  const CSCLayer* l2 = aState.aChamber->layer((*(aState.proto_segment.begin() + 1))->cscDetId().layer());
490  GlobalPoint gp2 = l2->toGlobal((*(aState.proto_segment.begin() + 1))->localPosition());
491  float err_stpos_h1 = (*(aState.proto_segment.begin()))->errorWithinStrip();
492  float err_stpos_h2 = (*(aState.proto_segment.begin() + 1))->errorWithinStrip();
493  const CSCLayer* l = aState.aChamber->layer(h->cscDetId().layer());
494  GlobalPoint hp = l->toGlobal(h->localPosition());
495  float err_stpos_h = h->errorWithinStrip();
496  float hphi = hp.phi(); // in (-pi, +pi]
497  if (hphi < 0.)
498  hphi += 2. * M_PI; // into range [0, 2pi)
499  float sphi = phiAtZ(aState, hp.z()); // in [0, 2*pi)
500  float phidif = sphi - hphi;
501  if (phidif < 0.)
502  phidif += 2. * M_PI; // into range [0, 2pi)
503  if (phidif > M_PI)
504  phidif -= 2. * M_PI; // into range (-pi, pi]
505  SVector6 r_glob;
506  CSCDetId id = h->cscDetId();
507  int iStn = id.iChamberType() - 1;
508  float dphi_incr = 0;
509  float pos_str = 1;
510  //increase dPhi cut if the hit is on the edge of the strip
511  float stpos = (*h).positionWithinStrip();
512  bool centr_str = false;
513  if (iStn != 0 && iStn != 1) {
514  if (stpos > -0.25 && stpos < 0.25)
515  centr_str = true;
516  }
517  if (err_stpos_h1 < 0.25 * strip_width[iStn] || err_stpos_h2 < 0.25 * strip_width[iStn] ||
518  err_stpos_h < 0.25 * strip_width[iStn]) {
519  dphi_incr = 0.5 * strip_width[iStn];
520  } else {
521  if (centr_str)
522  pos_str = 1.3;
523  }
524  r_glob((*(aState.proto_segment.begin()))->cscDetId().layer() - 1) = gp1.perp();
525  r_glob((*(aState.proto_segment.begin() + 1))->cscDetId().layer() - 1) = gp2.perp();
526  float R = hp.perp();
527  int layer = h->cscDetId().layer();
528  float r_interpolated = fit_r_phi(aState, r_glob, layer);
529  float dr = fabs(r_interpolated - R);
530  float maxWG_width[10] = {0, 0, 4.1, 5.69, 2.49, 5.06, 2.49, 5.06, 1.87, 5.06};
531  //find maxWG_width for ME11 (tilt = 29deg)
532  int wg_num = h->hitWire();
533  if (iStn == 0 || iStn == 1) {
534  if (wg_num == 1) {
535  maxWG_width[0] = 9.25;
536  maxWG_width[1] = 9.25;
537  }
538  if (wg_num > 1 && wg_num < 48) {
539  maxWG_width[0] = 3.14;
540  maxWG_width[1] = 3.14;
541  }
542  if (wg_num == 48) {
543  maxWG_width[0] = 10.75;
544  maxWG_width[1] = 10.75;
545  }
546  }
547 
548  if (aState.enlarge) {
549  return (fabs(phidif) < aState.dPhiIntMax * aState.strip_iadd * pos_str + dphi_incr &&
550  fabs(dr) < aState.dRIntMax * aState.strip_iadd * maxWG_width[iStn])
551  ? true
552  : false;
553 
554  } else {
555  return (fabs(phidif) < aState.dPhiIntMax * aState.strip_iadd * pos_str + dphi_incr &&
556  fabs(dr) < aState.dRIntMax * maxWG_width[iStn])
557  ? true
558  : false;
559  }
560 }

References CSCSegAlgoRU::AlgoState::aChamber, CSCSegAlgoRU::AlgoState::dPhiIntMax, flavorHistoryFilter_cfi::dr, CSCSegAlgoRU::AlgoState::dRIntMax, CSCSegAlgoRU::AlgoState::enlarge, fit_r_phi(), trackingPlots::hp, cmsLHEtoEOSManager::l, CSCChamber::layer(), M_PI, PV3DBase< T, PVType, FrameType >::perp(), phiAtZ(), CSCSegAlgoRU::AlgoState::proto_segment, dttmaxenums::R, CSCSegAlgoRU::AlgoState::strip_iadd, GeomDet::toGlobal(), and funct::true.

Referenced by tryAddingHitsToSegment().

◆ isSegmentGood()

bool CSCSegAlgoRU::isSegmentGood ( const AlgoState aState,
const ChamberHitContainer rechitsInChamber 
) const
private

Return true if segment is 'good'. In this algorithm, 'good' means has sufficient hits

Definition at line 577 of file CSCSegAlgoRU.cc.

577  {
578  // If the chamber has 20 hits or fewer, require at least 3 hits on segment
579  // If the chamber has >20 hits require at least 4 hits
580  // If it's the second cycle of the builder and there are <= 12 hits in chamber, require at least 3 hits on segment
581  //@@ THESE VALUES SHOULD BECOME PARAMETERS?
582  bool ok = false;
583  unsigned int iadd = (rechitsInChamber.size() > 20) ? 1 : 0;
584  if (aState.windowScale > 1.) {
585  iadd = 1;
586  if (rechitsInChamber.size() <= 12 && aState.enlarge)
587  iadd = 0;
588  }
589  if (aState.proto_segment.size() >= 3 + iadd)
590  ok = true;
591  return ok;
592 }

References CSCSegAlgoRU::AlgoState::enlarge, convertSQLiteXML::ok, CSCSegAlgoRU::AlgoState::proto_segment, and CSCSegAlgoRU::AlgoState::windowScale.

Referenced by buildSegments().

◆ phiAtZ()

float CSCSegAlgoRU::phiAtZ ( const AlgoState aState,
float  z 
) const
private

Always enforce direction of segment to point from IP outwards (Incorrect for particles not coming from IP, of course.)

Definition at line 562 of file CSCSegAlgoRU.cc.

562  {
563  if (!aState.sfit)
564  return 0.;
565  // Returns a phi in [ 0, 2*pi )
566  const CSCLayer* l1 = aState.aChamber->layer((*(aState.proto_segment.begin()))->cscDetId().layer());
567  GlobalPoint gp = l1->toGlobal(aState.sfit->intercept());
568  GlobalVector gv = l1->toGlobal(aState.sfit->localdir());
569  float x = gp.x() + (gv.x() / gv.z()) * (z - gp.z());
570  float y = gp.y() + (gv.y() / gv.z()) * (z - gp.z());
571  float phi = atan2(y, x);
572  if (phi < 0.f)
573  phi += 2. * M_PI;
574  return phi;
575 }

References CSCSegAlgoRU::AlgoState::aChamber, f, runTauDisplay::gp, CSCChamber::layer(), M_PI, phi, CSCSegAlgoRU::AlgoState::proto_segment, CSCSegAlgoRU::AlgoState::sfit, GeomDet::toGlobal(), x, PV3DBase< T, PVType, FrameType >::x(), y, PV3DBase< T, PVType, FrameType >::y(), z, and PV3DBase< T, PVType, FrameType >::z().

Referenced by isHitNearSegment().

◆ replaceHit()

bool CSCSegAlgoRU::replaceHit ( AlgoState aState,
const CSCRecHit2D h,
int  layer 
) const
private

Definition at line 887 of file CSCSegAlgoRU.cc.

887  {
888  // replace a hit from a layer
889  ChamberHitContainer::const_iterator it;
890  for (it = aState.proto_segment.begin(); it != aState.proto_segment.end();) {
891  if ((*it)->cscDetId().layer() == layer)
892  it = aState.proto_segment.erase(it);
893  else
894  ++it;
895  }
896  return addHit(aState, h, layer);
897 }

References addHit(), and CSCSegAlgoRU::AlgoState::proto_segment.

Referenced by compareProtoSegment().

◆ run()

std::vector<CSCSegment> CSCSegAlgoRU::run ( const CSCChamber aChamber,
const ChamberHitContainer rechits 
)
inlineoverride

Here we must implement the algorithm

Definition at line 77 of file CSCSegAlgoRU.h.

77  {
78  return buildSegments(aChamber, rechits);
79  }

References buildSegments(), and HI_PhotonSkim_cff::rechits.

◆ tryAddingHitsToSegment()

void CSCSegAlgoRU::tryAddingHitsToSegment ( AlgoState aState,
const ChamberHitContainer rechitsInChamber,
const BoolContainer used,
const LayerIndex layerIndex,
const ChamberHitContainerCIt  i1,
const ChamberHitContainerCIt  i2 
) const
private

Try adding non-used hits to segment

Definition at line 359 of file CSCSegAlgoRU.cc.

364  {
365  // Iterate over the layers with hits in the chamber
366  // Skip the layers containing the segment endpoints
367  // Test each hit on the other layers to see if it is near the segment
368  // If it is, see whether there is already a hit on the segment from the same layer
369  // - if so, and there are more than 2 hits on the segment, copy the segment,
370  // replace the old hit with the new hit. If the new segment chi2 is better
371  // then replace the original segment with the new one (by swap)
372  // - if not, copy the segment, add the hit. If the new chi2/dof is still satisfactory
373  // then replace the original segment with the new one (by swap)
375  ChamberHitContainerCIt ie = rechits.end();
376  for (ChamberHitContainerCIt i = ib; i != ie; ++i) {
377  int layer = layerIndex[i - ib];
378  if (hasHitOnLayer(aState, layer) && aState.proto_segment.size() <= 2)
379  continue;
380  if (layerIndex[i - ib] == layerIndex[i1 - ib] || layerIndex[i - ib] == layerIndex[i2 - ib] || used[i - ib])
381  continue;
382 
383  const CSCRecHit2D* h = *i;
384  if (isHitNearSegment(aState, h)) {
385  // Don't consider alternate hits on layers holding the two starting points
386  if (hasHitOnLayer(aState, layer)) {
387  if (aState.proto_segment.size() <= 2)
388  continue;
389  compareProtoSegment(aState, h, layer);
390  } else {
391  increaseProtoSegment(aState, h, layer, aState.chi2D_iadd);
392  }
393  } // h & seg close
394  } // i
395 }

References CSCSegAlgoRU::AlgoState::chi2D_iadd, compareProtoSegment(), hasHitOnLayer(), mps_fire::i, testProducerWithPsetDescEmpty_cfi::i1, testProducerWithPsetDescEmpty_cfi::i2, cuy::ib, increaseProtoSegment(), isHitNearSegment(), CSCSegAlgoRU::AlgoState::proto_segment, and HI_PhotonSkim_cff::rechits.

Referenced by buildSegments().

◆ updateParameters()

void CSCSegAlgoRU::updateParameters ( AlgoState aState) const
private

Definition at line 622 of file CSCSegAlgoRU.cc.

622  {
623  // Delete input CSCSegFit, create a new one and make the fit
624  // delete sfit;
625  aState.sfit = std::make_unique<CSCSegFit>(aState.aChamber, aState.proto_segment);
626  aState.sfit->fit();
627 }

References CSCSegAlgoRU::AlgoState::aChamber, CSCSegAlgoRU::AlgoState::proto_segment, and CSCSegAlgoRU::AlgoState::sfit.

Referenced by addHit(), and buildSegments().

Member Data Documentation

◆ chi2_str_

float CSCSegAlgoRU::chi2_str_
private

Definition at line 160 of file CSCSegAlgoRU.h.

Referenced by buildSegments(), and CSCSegAlgoRU().

◆ chi2Max

float CSCSegAlgoRU::chi2Max
private

Definition at line 159 of file CSCSegAlgoRU.h.

Referenced by buildSegments(), and CSCSegAlgoRU().

◆ chi2Norm_2D_

float CSCSegAlgoRU::chi2Norm_2D_
private

Definition at line 161 of file CSCSegAlgoRU.h.

Referenced by buildSegments(), and CSCSegAlgoRU().

◆ debugInfo

bool CSCSegAlgoRU::debugInfo
private

Definition at line 164 of file CSCSegAlgoRU.h.

◆ doCollisions

bool CSCSegAlgoRU::doCollisions
private

Definition at line 154 of file CSCSegAlgoRU.h.

Referenced by buildSegments(), and CSCSegAlgoRU().

◆ dPhiIntMax

float CSCSegAlgoRU::dPhiIntMax
private

Definition at line 158 of file CSCSegAlgoRU.h.

Referenced by buildSegments(), and CSCSegAlgoRU().

◆ dPhiMax

float CSCSegAlgoRU::dPhiMax
private

Definition at line 156 of file CSCSegAlgoRU.h.

Referenced by buildSegments(), and CSCSegAlgoRU().

◆ dRIntMax

float CSCSegAlgoRU::dRIntMax
private

Definition at line 157 of file CSCSegAlgoRU.h.

Referenced by buildSegments(), and CSCSegAlgoRU().

◆ dRMax

float CSCSegAlgoRU::dRMax
private

Definition at line 155 of file CSCSegAlgoRU.h.

Referenced by buildSegments(), and CSCSegAlgoRU().

◆ enlarge

bool CSCSegAlgoRU::enlarge
private

Definition at line 165 of file CSCSegAlgoRU.h.

Referenced by buildSegments(), and CSCSegAlgoRU().

◆ minLayersApart

int CSCSegAlgoRU::minLayersApart
private

Definition at line 163 of file CSCSegAlgoRU.h.

Referenced by CSCSegAlgoRU().

◆ myName

const std::string CSCSegAlgoRU::myName
private

Definition at line 148 of file CSCSegAlgoRU.h.

Referenced by CSCSegAlgoRU().

◆ theChi2

double CSCSegAlgoRU::theChi2
private

Definition at line 150 of file CSCSegAlgoRU.h.

◆ theDirection

LocalVector CSCSegAlgoRU::theDirection
private

Definition at line 152 of file CSCSegAlgoRU.h.

◆ theOrigin

LocalPoint CSCSegAlgoRU::theOrigin
private

Definition at line 151 of file CSCSegAlgoRU.h.

◆ uz

float CSCSegAlgoRU::uz
private

Definition at line 153 of file CSCSegAlgoRU.h.

◆ vz

float CSCSegAlgoRU::vz
private

Definition at line 153 of file CSCSegAlgoRU.h.

◆ wideSeg

float CSCSegAlgoRU::wideSeg
private

Definition at line 162 of file CSCSegAlgoRU.h.

Referenced by buildSegments(), and CSCSegAlgoRU().

Vector3DBase
Definition: Vector3DBase.h:8
GeomDet::position
const Surface::PositionType & position() const
The position (origin of the R.F.)
Definition: GeomDet.h:43
CSCSegAlgoRU::updateParameters
void updateParameters(AlgoState &aState) const
Definition: CSCSegAlgoRU.cc:622
DDAxes::y
testProducerWithPsetDescEmpty_cfi.i2
i2
Definition: testProducerWithPsetDescEmpty_cfi.py:46
CSCSegAlgoRU::replaceHit
bool replaceHit(AlgoState &aState, const CSCRecHit2D *h, int layer) const
Definition: CSCSegAlgoRU.cc:887
CSCSegAlgoRU::myName
const std::string myName
Definition: CSCSegAlgoRU.h:148
mps_fire.i
i
Definition: mps_fire.py:428
CSCChamber::layer
const CSCLayer * layer(CSCDetId id) const
Return the layer corresponding to the given id.
Definition: CSCChamber.cc:30
HLT_FULL_cff.points
points
Definition: HLT_FULL_cff.py:21521
f
double f[11][100]
Definition: MuScleFitUtils.cc:78
CSCSegAlgoRU::minLayersApart
int minLayersApart
Definition: CSCSegAlgoRU.h:163
PV3DBase::x
T x() const
Definition: PV3DBase.h:59
min
T min(T a, T b)
Definition: MathUtil.h:58
CSCSegAlgoRU::ChamberHitContainer
std::vector< const CSCRecHit2D * > ChamberHitContainer
Definition: CSCSegAlgoRU.h:51
CSCSegAlgoRU::tryAddingHitsToSegment
void tryAddingHitsToSegment(AlgoState &aState, const ChamberHitContainer &rechitsInChamber, const BoolContainer &used, const LayerIndex &layerIndex, const ChamberHitContainerCIt i1, const ChamberHitContainerCIt i2) const
Definition: CSCSegAlgoRU.cc:359
CSCDetId::ring
int ring() const
Definition: CSCDetId.h:68
testProducerWithPsetDescEmpty_cfi.i1
i1
Definition: testProducerWithPsetDescEmpty_cfi.py:45
CSCSegmentAlgorithmRU_cfi.chi2_str
chi2_str
Definition: CSCSegmentAlgorithmRU_cfi.py:7
CSCLayer
Definition: CSCLayer.h:24
CSCSegAlgoRU::chi2Max
float chi2Max
Definition: CSCSegAlgoRU.h:159
DDAxes::x
convertSQLiteXML.ok
bool ok
Definition: convertSQLiteXML.py:98
groupFilesInBlocks.temp
list temp
Definition: groupFilesInBlocks.py:142
CSCDetId::iChamberType
unsigned short iChamberType() const
Definition: CSCDetId.h:96
CSCSegAlgoRU::doCollisions
bool doCollisions
Definition: CSCSegAlgoRU.h:154
CSCSegAlgoRU::ChamberHitContainerCIt
std::vector< const CSCRecHit2D * >::const_iterator ChamberHitContainerCIt
Definition: CSCSegAlgoRU.h:52
CSCRecHit2D::cscDetId
CSCDetId cscDetId() const
Definition: CSCRecHit2D.h:58
groupFilesInBlocks.reverse
reverse
Definition: groupFilesInBlocks.py:131
edmScanValgrind.buffer
buffer
Definition: edmScanValgrind.py:171
testProducerWithPsetDescEmpty_cfi.z2
z2
Definition: testProducerWithPsetDescEmpty_cfi.py:41
CSCRecHit2D::errorWithinStrip
float errorWithinStrip() const
The uncertainty of the estimated position within the strip.
Definition: CSCRecHit2D.h:85
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
trackingPlots.hp
hp
Definition: trackingPlots.py:1248
SiPixelRawToDigiRegional_cfi.deltaPhi
deltaPhi
Definition: SiPixelRawToDigiRegional_cfi.py:9
CSCSegAlgoRU::baseline
void baseline(AlgoState &aState, int n_seg_min) const
Definition: CSCSegAlgoRU.cc:649
CSCSegAlgoRU::dPhiIntMax
float dPhiIntMax
Definition: CSCSegAlgoRU.h:158
CSCDetId::layer
int layer() const
Definition: CSCDetId.h:56
DDAxes::z
CSCSegAlgoRU::increaseProtoSegment
void increaseProtoSegment(AlgoState &aState, const CSCRecHit2D *h, int layer, int chi2_factor) const
Definition: CSCSegAlgoRU.cc:915
CSCSegAlgoRU::BoolContainer
std::vector< bool > BoolContainer
Definition: CSCSegAlgoRU.h:59
CSCSegAlgoRU::LayerIndex
std::vector< int > LayerIndex
Definition: CSCSegAlgoRU.h:50
HI_PhotonSkim_cff.rechits
rechits
Definition: HI_PhotonSkim_cff.py:76
h
dqmdumpme.k
k
Definition: dqmdumpme.py:60
Point3DBase< float, GlobalTag >
nhits
Definition: HIMultiTrackSelector.h:42
CSCSegment
Definition: CSCSegment.h:21
CSCSegAlgoRU::hasHitOnLayer
bool hasHitOnLayer(const AlgoState &aState, int layer) const
Definition: CSCSegAlgoRU.cc:878
runTauDisplay.gp
gp
Definition: runTauDisplay.py:431
funct::true
true
Definition: Factorize.h:173
GeomDet::toGlobal
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Definition: GeomDet.h:49
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:223
CSCRecHit2D
Definition: CSCRecHit2D.h:18
PV3DBase::barePhi
T barePhi() const
Definition: PV3DBase.h:65
CSCRecHit2D::localPosition
LocalPoint localPosition() const override
Definition: CSCRecHit2D.h:56
dumpMFGeometry_cfg.delta
delta
Definition: dumpMFGeometry_cfg.py:25
CSCDetId
Definition: CSCDetId.h:26
cand
Definition: decayParser.h:32
PV3DBase::y
T y() const
Definition: PV3DBase.h:60
CSCSegAlgoRU::flagHitsAsUsed
void flagHitsAsUsed(const AlgoState &aState, const ChamberHitContainer &rechitsInChamber, BoolContainer &used) const
Definition: CSCSegAlgoRU.cc:594
M_PI
#define M_PI
Definition: BXVectorInputProducer.cc:49
CSCSegAlgoRU::enlarge
bool enlarge
Definition: CSCSegAlgoRU.h:165
cuy.ib
ib
Definition: cuy.py:662
cmsLHEtoEOSManager.l
l
Definition: cmsLHEtoEOSManager.py:204
CSCSegAlgoRU::isSegmentGood
bool isSegmentGood(const AlgoState &aState, const ChamberHitContainer &rechitsInChamber) const
Definition: CSCSegAlgoRU.cc:577
DDAxes::phi
hi
Definition: EPCuts.h:4
CSCSegAlgoRU::areHitsCloseInGlobalPhi
bool areHitsCloseInGlobalPhi(const AlgoState &aState, const CSCRecHit2D *h1, const CSCRecHit2D *h2) const
Definition: CSCSegAlgoRU.cc:444
CSCRecHit2D::channels
int channels(unsigned int i) const
Extracting strip channel numbers comprising the rechit - low.
Definition: CSCRecHit2D.h:61
eostools.move
def move(src, dest)
Definition: eostools.py:511
CSCSegAlgoRU::dRIntMax
float dRIntMax
Definition: CSCSegAlgoRU.h:157
CSCSegAlgoRU::dPhiMax
float dPhiMax
Definition: CSCSegAlgoRU.h:156
flavorHistoryFilter_cfi.dr
dr
Definition: flavorHistoryFilter_cfi.py:37
CSCSegAlgoRU::fitX
float fitX(const AlgoState &aState, SVector6 points, SVector6 errors, int ir, int ir2, float &chi2_str) const
Definition: CSCSegAlgoRU.cc:844
CSCSegAlgoRU::addHit
bool addHit(AlgoState &aState, const CSCRecHit2D *hit, int layer) const
Utility functions.
Definition: CSCSegAlgoRU.cc:608
CSCSegAlgoRU::buildSegments
std::vector< CSCSegment > buildSegments(const CSCChamber *aChamber, const ChamberHitContainer &rechits) const
Definition: CSCSegAlgoRU.cc:57
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
S
Definition: CSCDBL1TPParametersExtended.h:16
CSCSegmentAlgorithm::CSCSegmentAlgorithm
CSCSegmentAlgorithm(const edm::ParameterSet &)
Constructor.
Definition: CSCSegmentAlgorithm.h:26
CSCSegAlgoRU::chi2Norm_2D_
float chi2Norm_2D_
Definition: CSCSegAlgoRU.h:161
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
CSCSegAlgoRU::SVector6
ROOT::Math::SVector< double, 6 > SVector6
Typedefs.
Definition: CSCSegAlgoRU.h:48
CSCSegAlgoRU::dRMax
float dRMax
Definition: CSCSegAlgoRU.h:155
CSCDetId::station
int station() const
Definition: CSCDetId.h:79
CSCRecHit2D::nStrips
unsigned int nStrips() const
Definition: CSCRecHit2D.h:62
CSCRecHit2D::hitWire
short int hitWire() const
L1A.
Definition: CSCRecHit2D.h:68
dqmiolumiharvest.j
j
Definition: dqmiolumiharvest.py:66
slope
static const double slope[3]
Definition: CastorTimeSlew.cc:6
CSCSegAlgoRU::phiAtZ
float phiAtZ(const AlgoState &aState, float z) const
Definition: CSCSegAlgoRU.cc:562
dttmaxenums::R
Definition: DTTMax.h:29
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:69
CSCSegAlgoRU::isHitNearSegment
bool isHitNearSegment(const AlgoState &aState, const CSCRecHit2D *h) const
Definition: CSCSegAlgoRU.cc:472
CSCSegAlgoRU::fit_r_phi
float fit_r_phi(const AlgoState &aState, const SVector6 &points, int layer) const
Definition: CSCSegAlgoRU.cc:629
CSCSegAlgoRU::areHitsCloseInR
bool areHitsCloseInR(const AlgoState &aState, const CSCRecHit2D *h1, const CSCRecHit2D *h2) const
Utility functions.
Definition: CSCSegAlgoRU.cc:397
kLayer
static const std::string kLayer("layer")
CSCSegAlgoRU::compareProtoSegment
void compareProtoSegment(AlgoState &aState, const CSCRecHit2D *h, int layer) const
Definition: CSCSegAlgoRU.cc:899
CSCSegAlgoRU::chi2_str_
float chi2_str_
Definition: CSCSegAlgoRU.h:160
CSCSegAlgoRU::wideSeg
float wideSeg
Definition: CSCSegAlgoRU.h:162
debug_messages_cfi.errors
errors
Definition: debug_messages_cfi.py:54