00001 #include "RecoTracker/TkDetLayers/interface/PixelBlade.h"
00002
00003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00004
00005 #include "RecoTracker/TkDetLayers/interface/BladeShapeBuilderFromDet.h"
00006 #include "RecoTracker/TkDetLayers/interface/LayerCrossingSide.h"
00007 #include "RecoTracker/TkDetLayers/interface/DetGroupMerger.h"
00008 #include "RecoTracker/TkDetLayers/interface/CompatibleDetToGroupAdder.h"
00009
00010 #include "TrackingTools/DetLayers/interface/DetLayerException.h"
00011 #include "TrackingTools/PatternTools/interface/MeasurementEstimator.h"
00012 #include "TrackingTools/GeomPropagators/interface/HelixArbitraryPlaneCrossing.h"
00013
00014 using namespace std;
00015
00016 typedef GeometricSearchDet::DetWithState DetWithState;
00017
00018 PixelBlade::PixelBlade(vector<const GeomDet*>& frontDets,
00019 vector<const GeomDet*>& backDets):
00020 theFrontDets(frontDets), theBackDets(backDets)
00021 {
00022 theDets.assign(theFrontDets.begin(),theFrontDets.end());
00023 theDets.insert(theDets.end(),theBackDets.begin(),theBackDets.end());
00024
00025 theDiskSector = BladeShapeBuilderFromDet()(theDets);
00026 theFrontDiskSector = BladeShapeBuilderFromDet()(theFrontDets);
00027 theBackDiskSector = BladeShapeBuilderFromDet()(theBackDets);
00028
00029
00030
00031 LogDebug("TkDetLayers") << "DEBUG INFO for PixelBlade" ;
00032 LogDebug("TkDetLayers") << "Blade z, perp, innerRadius, outerR: "
00033 << this->position().z() << " , "
00034 << this->position().perp() << " , "
00035 << theDiskSector->innerRadius() << " , "
00036 << theDiskSector->outerRadius() ;
00037
00038 for(vector<const GeomDet*>::const_iterator it=theFrontDets.begin();
00039 it!=theFrontDets.end(); it++){
00040 LogDebug("TkDetLayers") << "frontDet phi,z,r: "
00041 << (*it)->position().phi() << " , "
00042 << (*it)->position().z() << " , "
00043 << (*it)->position().perp() ;;
00044 }
00045
00046 for(vector<const GeomDet*>::const_iterator it=theBackDets.begin();
00047 it!=theBackDets.end(); it++){
00048 LogDebug("TkDetLayers") << "backDet phi,z,r: "
00049 << (*it)->position().phi() << " , "
00050 << (*it)->position().z() << " , "
00051 << (*it)->position().perp() ;
00052 }
00053
00054
00055 }
00056
00057
00058 const vector<const GeometricSearchDet*>&
00059 PixelBlade::components() const{
00060 throw DetLayerException("TOBRod doesn't have GeometricSearchDet components");
00061 }
00062
00063 pair<bool, TrajectoryStateOnSurface>
00064 PixelBlade::compatible( const TrajectoryStateOnSurface& ts, const Propagator&,
00065 const MeasurementEstimator&) const{
00066 edm::LogError("TkDetLayers") << "temporary dummy implementation of PixelBlade::compatible()!!" ;
00067 return pair<bool,TrajectoryStateOnSurface>();
00068 }
00069
00070
00071
00072 void
00073 PixelBlade::groupedCompatibleDetsV( const TrajectoryStateOnSurface& tsos,
00074 const Propagator& prop,
00075 const MeasurementEstimator& est,
00076 std::vector<DetGroup> & result) const{
00077 SubLayerCrossings crossings;
00078 crossings = computeCrossings( tsos, prop.propagationDirection());
00079 if(! crossings.isValid()) return;
00080
00081 vector<DetGroup> closestResult;
00082 addClosest( tsos, prop, est, crossings.closest(), closestResult);
00083
00084 if (closestResult.empty()){
00085 vector<DetGroup> nextResult;
00086 addClosest( tsos, prop, est, crossings.other(), nextResult);
00087 if(nextResult.empty()) return;
00088
00089 DetGroupElement nextGel( nextResult.front().front());
00090 int crossingSide = LayerCrossingSide().endcapSide( nextGel.trajectoryState(), prop);
00091
00092 DetGroupMerger::orderAndMergeTwoLevels( closestResult, nextResult, result,
00093 crossings.closestIndex(), crossingSide);
00094 }
00095 else {
00096 DetGroupElement closestGel( closestResult.front().front());
00097 float window = computeWindowSize( closestGel.det(), closestGel.trajectoryState(), est);
00098
00099 searchNeighbors( tsos, prop, est, crossings.closest(), window,
00100 closestResult, false);
00101
00102 vector<DetGroup> nextResult;
00103 searchNeighbors( tsos, prop, est, crossings.other(), window,
00104 nextResult, true);
00105
00106 int crossingSide = LayerCrossingSide().endcapSide( closestGel.trajectoryState(), prop);
00107 DetGroupMerger::orderAndMergeTwoLevels( closestResult, nextResult, result,
00108 crossings.closestIndex(), crossingSide);
00109 }
00110 }
00111
00112 SubLayerCrossings
00113 PixelBlade::computeCrossings( const TrajectoryStateOnSurface& startingState,
00114 PropagationDirection propDir) const
00115 {
00116 HelixPlaneCrossing::PositionType startPos( startingState.globalPosition());
00117 HelixPlaneCrossing::DirectionType startDir( startingState.globalMomentum());
00118 double rho( startingState.transverseCurvature());
00119
00120 HelixArbitraryPlaneCrossing crossing( startPos, startDir, rho, propDir);
00121
00122 pair<bool,double> innerPath = crossing.pathLength( *theFrontDiskSector);
00123 if (!innerPath.first) return SubLayerCrossings();
00124
00125 GlobalPoint gInnerPoint( crossing.position(innerPath.second));
00126
00127
00128
00129 int innerIndex = findBin(gInnerPoint.perp(),0);
00130 float innerDist = fabs( findPosition(innerIndex,0).perp() - gInnerPoint.perp());
00131 SubLayerCrossing innerSLC( 0, innerIndex, gInnerPoint);
00132
00133 pair<bool,double> outerPath = crossing.pathLength( *theBackDiskSector);
00134 if (!outerPath.first) return SubLayerCrossings();
00135
00136 GlobalPoint gOuterPoint( crossing.position(outerPath.second));
00137
00138
00139
00140 int outerIndex = findBin(gOuterPoint.perp(),1);
00141 float outerDist = fabs( findPosition(outerIndex,1).perp() - gOuterPoint.perp());
00142 SubLayerCrossing outerSLC( 1, outerIndex, gOuterPoint);
00143
00144 if (innerDist < outerDist) {
00145 return SubLayerCrossings( innerSLC, outerSLC, 0);
00146 }
00147 else {
00148 return SubLayerCrossings( outerSLC, innerSLC, 1);
00149 }
00150 }
00151
00152
00153
00154
00155 bool
00156 PixelBlade::addClosest( const TrajectoryStateOnSurface& tsos,
00157 const Propagator& prop,
00158 const MeasurementEstimator& est,
00159 const SubLayerCrossing& crossing,
00160 vector<DetGroup>& result) const
00161 {
00162
00163 const vector<const GeomDet*>& sBlade( subBlade( crossing.subLayerIndex()));
00164 return CompatibleDetToGroupAdder().add( *sBlade[crossing.closestDetIndex()],
00165 tsos, prop, est, result);
00166 }
00167
00168
00169 float PixelBlade::computeWindowSize( const GeomDet* det,
00170 const TrajectoryStateOnSurface& tsos,
00171 const MeasurementEstimator& est) const
00172 {
00173 return
00174 est.maximalLocalDisplacement(tsos, det->surface()).x();
00175 }
00176
00177
00178
00179
00180 void PixelBlade::searchNeighbors( const TrajectoryStateOnSurface& tsos,
00181 const Propagator& prop,
00182 const MeasurementEstimator& est,
00183 const SubLayerCrossing& crossing,
00184 float window,
00185 vector<DetGroup>& result,
00186 bool checkClosest) const
00187 {
00188 GlobalPoint gCrossingPos = crossing.position();
00189
00190 const vector<const GeomDet*>& sBlade( subBlade( crossing.subLayerIndex()));
00191
00192 int closestIndex = crossing.closestDetIndex();
00193 int negStartIndex = closestIndex-1;
00194 int posStartIndex = closestIndex+1;
00195
00196 if (checkClosest) {
00197 if (gCrossingPos.perp() < sBlade[closestIndex]->surface().position().perp()) {
00198 posStartIndex = closestIndex;
00199 }
00200 else {
00201 negStartIndex = closestIndex;
00202 }
00203 }
00204
00205 typedef CompatibleDetToGroupAdder Adder;
00206 for (int idet=negStartIndex; idet >= 0; idet--) {
00207 if (!overlap( gCrossingPos, *sBlade[idet], window)) break;
00208 if (!Adder::add( *sBlade[idet], tsos, prop, est, result)) break;
00209 }
00210 for (int idet=posStartIndex; idet < static_cast<int>(sBlade.size()); idet++) {
00211 if (!overlap( gCrossingPos, *sBlade[idet], window)) break;
00212 if (!Adder::add( *sBlade[idet], tsos, prop, est, result)) break;
00213 }
00214 }
00215
00216
00217
00218 bool PixelBlade::overlap( const GlobalPoint& crossPoint, const GeomDet& det, float window) const
00219 {
00220
00221
00222
00223 const float relativeMargin = 1.01;
00224
00225 LocalPoint localCrossPoint( det.surface().toLocal(crossPoint));
00226
00227
00228
00229
00230
00231 float localX = localCrossPoint.x();
00232 float detHalfLength = det.surface().bounds().length()/2.;
00233
00234
00235
00236
00237 if ( ( fabs(localX)-window) < relativeMargin*detHalfLength ) {
00238 return true;
00239 } else {
00240 return false;
00241 }
00242 }
00243
00244 int
00245 PixelBlade::findBin( float R,int diskSectorIndex) const
00246 {
00247 vector<const GeomDet*> localDets = diskSectorIndex==0 ? theFrontDets : theBackDets;
00248
00249 int theBin = 0;
00250 float rDiff = fabs( R - localDets.front()->surface().position().perp());;
00251 for (vector<const GeomDet*>::const_iterator i=localDets.begin(); i !=localDets.end(); i++){
00252 float testDiff = fabs( R - (**i).surface().position().perp());
00253 if ( testDiff < rDiff) {
00254 rDiff = testDiff;
00255 theBin = i - localDets.begin();
00256 }
00257 }
00258 return theBin;
00259 }
00260
00261
00262
00263 GlobalPoint
00264 PixelBlade::findPosition(int index,int diskSectorType) const
00265 {
00266 vector<const GeomDet*> diskSector = diskSectorType == 0 ? theFrontDets : theBackDets;
00267 return (diskSector[index])->surface().position();
00268 }
00269