TkDetLayers/src/Phase2EndcapLayerDoubleDisk.cc

0001 #include "Phase2EndcapLayerDoubleDisk.h"
0002
0003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0004
0005 #include "DataFormats/GeometrySurface/interface/SimpleDiskBounds.h"
0006
0007 #include "TrackingTools/DetLayers/interface/DetLayerException.h"
0008 #include "TrackingTools/GeomPropagators/interface/HelixForwardPlaneCrossing.h"
0009
0010 #include <array>
0011 #include "DetGroupMerger.h"
0012
0013 using namespace std;
0014
0015 const std::vector<const GeometricSearchDet*>& Phase2EndcapLayerDoubleDisk::components() const {
0016   throw cms::Exception("Phase2EndcapLayerDoubleDisk::components() is not implemented");
0017 }
0018
0019 void Phase2EndcapLayerDoubleDisk::fillSubDiskPars(int i) {
0020   const BoundDisk& subDiskDisk = static_cast<const BoundDisk&>(theComps[i]->surface());
0021   SubDiskPar tempPar;
0022   tempPar.theSubDiskZ = std::abs(subDiskDisk.position().z());
0023   subDiskPars.push_back(tempPar);
0024 }
0025
0026 Phase2EndcapLayerDoubleDisk::Phase2EndcapLayerDoubleDisk(vector<const Phase2EndcapSubDisk*>& subDisks)
0027     : RingedForwardLayer(true), theComponents{nullptr} {
0028   theSubDisksSize = subDisks.size();
0029   LogDebug("TkDetLayers") << "Number of subdisks in Phase2 IT EC layer is " << theSubDisksSize << std::endl;
0030   setSurface(computeDisk(subDisks));
0031
0032   for (unsigned int i = 0; i != subDisks.size(); ++i) {
0033     theComps.push_back(subDisks[i]);
0034     fillSubDiskPars(i);
0035     theBasicComps.insert(
0036         theBasicComps.end(), (*subDisks[i]).basicComponents().begin(), (*subDisks[i]).basicComponents().end());
0037   }
0038
0039   LogDebug("TkDetLayers") << "==== DEBUG Phase2EndcapLayer =====";
0040   LogDebug("TkDetLayers") << "r,zed pos  , thickness, innerR, outerR: " << this->position().perp() << " , "
0041                           << this->position().z() << " , " << this->specificSurface().bounds().thickness() << " , "
0042                           << this->specificSurface().innerRadius() << " , " << this->specificSurface().outerRadius();
0043 }
0044
0045 BoundDisk* Phase2EndcapLayerDoubleDisk::computeDisk(const vector<const Phase2EndcapSubDisk*>& subDisks) const {
0046   return tkDetUtil::computeDisk(subDisks);
0047 }
0048
0049 Phase2EndcapLayerDoubleDisk::~Phase2EndcapLayerDoubleDisk() {
0050   for (auto c : theComps)
0051     delete c;
0052
0053   delete theComponents.load();
0054 }
0055
0056 void Phase2EndcapLayerDoubleDisk::groupedCompatibleDetsV(const TrajectoryStateOnSurface& startingState,
0057                                                          const Propagator& prop,
0058                                                          const MeasurementEstimator& est,
0059                                                          std::vector<DetGroup>& result) const {
0060   std::array<int, 2> const& subDiskIndices = subDiskIndicesByCrossingProximity(startingState, prop);
0061
0062   //order subdisks in z
0063   //Subdisk near in z: 0, Subdisk far in z: 1
0064   std::vector<int> subDiskOrder(theSubDisksSize);
0065   std::fill(subDiskOrder.begin(), subDiskOrder.end(), 1);
0066   if (theSubDisksSize > 1) {
0067     if (std::abs(theComps[0]->position().z()) < std::abs(theComps[1]->position().z())) {
0068       for (int i = 0; i < theSubDisksSize; i++) {
0069         if (i % 2 == 0)
0070           subDiskOrder[i] = 0;
0071       }
0072     } else if (std::abs(theComps[0]->position().z()) > std::abs(theComps[1]->position().z())) {
0073       std::fill(subDiskOrder.begin(), subDiskOrder.end(), 0);
0074       for (int i = 0; i < theSubDisksSize; i++) {
0075         if (i % 2 == 0)
0076           subDiskOrder[i] = 1;
0077       }
0078     } else {
0079       throw DetLayerException("SubDisks in Endcap Layer have same z position, no idea how to order them!");
0080     }
0081   }
0082
0083   auto index = [&subDiskIndices, &subDiskOrder](int i) { return subDiskOrder[subDiskIndices[i]]; };
0084
0085   std::vector<DetGroup> closestResult;
0086   theComps[subDiskIndices[0]]->groupedCompatibleDetsV(startingState, prop, est, closestResult);
0087   // if the closest is empty, use the next one and exit: inherited from TID !
0088   if (closestResult.empty()) {
0089     theComps[subDiskIndices[1]]->groupedCompatibleDetsV(startingState, prop, est, result);
0090     return;
0091   }
0092
0093   // check if next subdisk is found
0094
0095   bool subdisk1ok = subDiskIndices[1] != -1;
0096
0097   // determine if we are propagating from in to out (0) or from out to in (1)
0098
0099   int direction = 0;
0100   if (startingState.globalPosition().z() * startingState.globalMomentum().z() > 0) {
0101     if (prop.propagationDirection() == alongMomentum)
0102       direction = 0;
0103     else
0104       direction = 1;
0105   } else {
0106     if (prop.propagationDirection() == alongMomentum)
0107       direction = 1;
0108     else
0109       direction = 0;
0110   }
0111
0112   if (index(0) == index(1)) {
0113     if (subdisk1ok) {
0114       std::vector<DetGroup> subdisk1res;
0115       theComps[subDiskIndices[1]]->groupedCompatibleDetsV(startingState, prop, est, subdisk1res);
0116       DetGroupMerger::addSameLevel(std::move(subdisk1res), closestResult);
0117       result.swap(closestResult);
0118       return;
0119     }
0120   } else {
0121     std::vector<DetGroup> subdisk1res;
0122     if (subdisk1ok) {
0123       theComps[subDiskIndices[1]]->groupedCompatibleDetsV(startingState, prop, est, subdisk1res);
0124     }
0125     if (!subdisk1res.empty()) {
0126       DetGroupMerger::orderAndMergeTwoLevels(
0127           std::move(closestResult), std::move(subdisk1res), result, index(0), direction);
0128       return;
0129     } else {
0130       result.swap(closestResult);
0131       return;
0132     }
0133   }
0134 }
0135
0136 std::array<int, 2> Phase2EndcapLayerDoubleDisk::subDiskIndicesByCrossingProximity(
0137     const TrajectoryStateOnSurface& startingState, const Propagator& prop) const {
0138   typedef HelixForwardPlaneCrossing Crossing;
0139   typedef MeasurementEstimator::Local2DVector Local2DVector;
0140
0141   HelixPlaneCrossing::PositionType startPos(startingState.globalPosition());
0142   HelixPlaneCrossing::DirectionType startDir(startingState.globalMomentum());
0143   PropagationDirection propDir(prop.propagationDirection());
0144   float rho(startingState.transverseCurvature());
0145
0146   // calculate the crossings with the subdisk surfaces
0147
0148   Crossing myXing(startPos, startDir, rho, propDir);
0149
0150   std::vector<GlobalPoint> subDiskCrossings;
0151   subDiskCrossings.reserve(theSubDisksSize);
0152
0153   for (int i = 0; i < theSubDisksSize; i++) {
0154     const BoundDisk& theSubDisk = static_cast<const BoundDisk&>(theComps[i]->surface());
0155     pair<bool, double> pathlen = myXing.pathLength(theSubDisk);
0156     if (pathlen.first) {
0157       subDiskCrossings.push_back(GlobalPoint(myXing.position(pathlen.second)));
0158     } else {
0159       // TO FIX.... perhaps there is something smarter to do
0160       subDiskCrossings.push_back(GlobalPoint(0., 0., 0.));
0161     }
0162   }
0163
0164   //find two closest subdisks to the crossing
0165
0166   return findTwoClosest(subDiskCrossings);
0167 }
0168
0169 std::array<int, 2> Phase2EndcapLayerDoubleDisk::findTwoClosest(std::vector<GlobalPoint> subDiskCrossing) const {
0170   std::array<int, 2> theBins = {{-1, -1}};
0171   theBins[0] = 0;
0172   float initialZ = subDiskPars[0].theSubDiskZ;
0173   float zDiff0 = std::abs(subDiskCrossing[0].z() - initialZ);
0174   float zDiff1 = -1.;
0175   for (int i = 1; i < theSubDisksSize; i++) {
0176     float subDiskZ = subDiskPars[i].theSubDiskZ;
0177     float testDiff = std::abs(subDiskCrossing[i].z() - subDiskZ);
0178     if (testDiff < zDiff0) {
0179       zDiff1 = zDiff0;
0180       zDiff0 = testDiff;
0181       theBins[1] = theBins[0];
0182       theBins[0] = i;
0183     } else if (zDiff1 < 0 || testDiff < zDiff1) {
0184       zDiff1 = testDiff;
0185       theBins[1] = i;
0186     }
0187   }
0188
0189   return theBins;
0190 }