00001 #include "RecoTracker/SpecialSeedGenerators/interface/CtfSpecialSeedGenerator.h"
00002
00003 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
00004 #include "Geometry/Records/interface/IdealGeometryRecord.h"
00005 #include "DataFormats/GeometrySurface/interface/RectangularPlaneBounds.h"
00006 #include "TrackingTools/GeomPropagators/interface/StraightLinePlaneCrossing.h"
00007
00008 #include "RecoTracker/TkTrackingRegions/interface/TrackingRegionProducerFactory.h"
00009 #include "RecoTracker/TkTrackingRegions/interface/TrackingRegion.h"
00010 #include "RecoTracker/TkTrackingRegions/interface/OrderedHitsGeneratorFactory.h"
00011 #include "RecoTracker/TkSeedingLayers/interface/OrderedSeedingHits.h"
00012
00013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00014
00015 #include "RecoTracker/SpecialSeedGenerators/interface/ClusterChecker.h"
00016
00017 using namespace ctfseeding;
00018
00019 CtfSpecialSeedGenerator::CtfSpecialSeedGenerator(const edm::ParameterSet& conf):
00020 conf_(conf),
00021 requireBOFF(conf.getParameter<bool>("requireBOFF")),
00022 theMaxSeeds(conf.getParameter<int32_t>("maxSeeds"))
00023 {
00024 useScintillatorsConstraint = conf_.getParameter<bool>("UseScintillatorsConstraint");
00025 edm::LogVerbatim("CtfSpecialSeedGenerator") << "Constructing CtfSpecialSeedGenerator";
00026 produces<TrajectorySeedCollection>();
00027 theSeedBuilder =0;
00028 theRegionProducer =0;
00029 }
00030
00031 CtfSpecialSeedGenerator::~CtfSpecialSeedGenerator(){
00032 }
00033
00034 void CtfSpecialSeedGenerator::endRun(edm::Run const&, edm::EventSetup const&){
00035 if (theSeedBuilder) { delete theSeedBuilder; theSeedBuilder = 0; }
00036 if (theRegionProducer) { delete theRegionProducer; theRegionProducer = 0; }
00037 std::vector<OrderedHitsGenerator*>::iterator iGen;
00038 for (iGen = theGenerators.begin(); iGen != theGenerators.end(); iGen++){
00039 delete (*iGen);
00040 }
00041 theGenerators.clear();
00042 }
00043
00044 void CtfSpecialSeedGenerator::beginRun(edm::Run const&, const edm::EventSetup& iSetup){
00045 std::string builderName = conf_.getParameter<std::string>("TTRHBuilder");
00046 iSetup.get<TransientRecHitRecord>().get(builderName,theBuilder);
00047
00048 iSetup.get<IdealMagneticFieldRecord>().get(theMagfield);
00049 iSetup.get<TrackerDigiGeometryRecord>().get(theTracker);
00050
00051 edm::LogVerbatim("CtfSpecialSeedGenerator") << "Initializing...";
00052 if (useScintillatorsConstraint){
00053 edm::ParameterSet upperScintPar = conf_.getParameter<edm::ParameterSet>("UpperScintillatorParameters");
00054 edm::ParameterSet lowerScintPar = conf_.getParameter<edm::ParameterSet>("LowerScintillatorParameters");
00055 RectangularPlaneBounds upperBounds(upperScintPar.getParameter<double>("WidthInX"),
00056 upperScintPar.getParameter<double>("LenghtInZ"),
00057 1);
00058 GlobalPoint upperPosition(upperScintPar.getParameter<double>("GlobalX"),
00059 upperScintPar.getParameter<double>("GlobalY"),
00060 upperScintPar.getParameter<double>("GlobalZ"));
00061 edm::LogVerbatim("CtfSpecialSeedGenerator")
00062 << "Upper Scintillator position x, y, z " << upperPosition.x()
00063 << ", " << upperPosition.y() << ", " << upperPosition.z();
00064 RectangularPlaneBounds lowerBounds(lowerScintPar.getParameter<double>("WidthInX"),
00065 lowerScintPar.getParameter<double>("LenghtInZ"),
00066 1);
00067 GlobalPoint lowerPosition(lowerScintPar.getParameter<double>("GlobalX"),
00068 lowerScintPar.getParameter<double>("GlobalY"),
00069 lowerScintPar.getParameter<double>("GlobalZ"));
00070 edm::LogVerbatim("CtfSpecialSeedGenerator")
00071 << "Lower Scintillator position x, y, z " << lowerPosition.x()
00072 << ", " << lowerPosition.y() << ", " << lowerPosition.z() ;
00073 TkRotation<float> rot(1,0,0,0,0,1,0,1,0);
00074 upperScintillator = BoundPlane::build(upperPosition, rot, &upperBounds);
00075 lowerScintillator = BoundPlane::build(lowerPosition, rot, &lowerBounds);
00076 }
00077 edm::ParameterSet regfactoryPSet = conf_.getParameter<edm::ParameterSet>("RegionFactoryPSet");
00078 std::string regfactoryName = regfactoryPSet.getParameter<std::string>("ComponentName");
00079 theRegionProducer = TrackingRegionProducerFactory::get()->create(regfactoryName,regfactoryPSet);
00080
00081 edm::ESHandle<Propagator> propagatorAlongHandle;
00082 iSetup.get<TrackingComponentsRecord>().get("PropagatorWithMaterial",propagatorAlongHandle);
00083 edm::ESHandle<Propagator> propagatorOppositeHandle;
00084 iSetup.get<TrackingComponentsRecord>().get("PropagatorWithMaterialOpposite",propagatorOppositeHandle);
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107 std::vector<edm::ParameterSet> pSets = conf_.getParameter<std::vector<edm::ParameterSet> >("OrderedHitsFactoryPSets");
00108 std::vector<edm::ParameterSet>::const_iterator iPSet;
00109 for (iPSet = pSets.begin(); iPSet != pSets.end(); iPSet++){
00110 std::string hitsfactoryName = iPSet->getParameter<std::string>("ComponentName");
00111 theGenerators.push_back(OrderedHitsGeneratorFactory::get()->create( hitsfactoryName, *iPSet));
00112 std::string propagationDirection = iPSet->getParameter<std::string>("PropagationDirection");
00113 if (propagationDirection == "alongMomentum") thePropDirs.push_back(alongMomentum);
00114 else thePropDirs.push_back(oppositeToMomentum);
00115 std::string navigationDirection = iPSet->getParameter<std::string>("NavigationDirection");
00116 if (navigationDirection == "insideOut") theNavDirs.push_back(insideOut);
00117 else theNavDirs.push_back(outsideIn);
00118 edm::LogVerbatim("CtfSpecialSeedGenerator") << "hitsGenerator done";
00119 }
00120 bool setMomentum = conf_.getParameter<bool>("SetMomentum");
00121 std::vector<int> charges;
00122 if (setMomentum){
00123 charges = conf_.getParameter<std::vector<int> >("Charges");
00124 }
00125 theSeedBuilder = new SeedFromGenericPairOrTriplet(theMagfield.product(),
00126 theTracker.product(),
00127 theBuilder.product(),
00128 propagatorAlongHandle.product(),
00129 propagatorOppositeHandle.product(),
00130 charges,
00131 setMomentum,
00132 conf_.getParameter<double>("ErrorRescaling"));
00133 double p = 1;
00134 if (setMomentum) {
00135 p = conf_.getParameter<double>("SeedMomentum");
00136 theSeedBuilder->setMomentumTo(p);
00137 }
00138
00139 }
00140
00141 void CtfSpecialSeedGenerator::produce(edm::Event& e, const edm::EventSetup& iSetup)
00142 {
00143
00144 std::auto_ptr<TrajectorySeedCollection> output(new TrajectorySeedCollection);
00145
00146
00147 ClusterChecker check(conf_);
00148 if ( !requireBOFF || (theMagfield->inTesla(GlobalPoint(0,0,0)).mag() == 0.00) ) {
00149 size_t clustsOrZero = check.tooManyClusters(e);
00150 if (!clustsOrZero){
00151 bool ok = run(iSetup, e, *output);
00152 if (!ok) { ; }
00153 } else edm::LogError("TooManyClusters") << "Found too many clusters (" << clustsOrZero << "), bailing out.\n";
00154 }
00155
00156
00157 edm::LogVerbatim("CtfSpecialSeedGenerator") << " number of seeds = "<< output->size();
00158 e.put(output);
00159 }
00160
00161 bool CtfSpecialSeedGenerator::run(const edm::EventSetup& iSetup,
00162 const edm::Event& e,
00163 TrajectorySeedCollection& output){
00164 std::vector<TrackingRegion*> regions = theRegionProducer->regions(e, iSetup);
00165 std::vector<TrackingRegion*>::const_iterator iReg;
00166 bool ok = true;
00167 for (iReg = regions.begin(); iReg != regions.end(); iReg++){
00168 if(!theSeedBuilder->momentumFromPSet()) theSeedBuilder->setMomentumTo((*iReg)->ptMin());
00169 std::vector<OrderedHitsGenerator*>::const_iterator iGen;
00170 int i = 0;
00171 for (iGen = theGenerators.begin(); iGen != theGenerators.end(); iGen++){
00172 ok = buildSeeds(iSetup,
00173 e,
00174 (*iGen)->run(**iReg, e, iSetup),
00175 theNavDirs[i],
00176 thePropDirs[i],
00177 output);
00178 i++;
00179 if (!ok) break;
00180 }
00181 if (!ok) break;
00182 }
00183
00184 for (std::vector<TrackingRegion*>::iterator iReg = regions.begin(); iReg != regions.end(); iReg++){
00185 delete *iReg;
00186 }
00187 return ok;
00188 }
00189
00190 bool CtfSpecialSeedGenerator::buildSeeds(const edm::EventSetup& iSetup,
00191 const edm::Event& e,
00192 const OrderedSeedingHits& osh,
00193 const NavigationDirection& navdir,
00194 const PropagationDirection& dir,
00195 TrajectorySeedCollection& output){
00196
00197 edm::LogInfo("CtfSpecialSeedGenerator")<<"osh.size() " << osh.size();
00198 for (unsigned int i = 0; i < osh.size(); i++){
00199 SeedingHitSet shs = osh[i];
00200 if (preliminaryCheck(shs,iSetup)){
00201 std::vector<TrajectorySeed*> seeds = theSeedBuilder->seed(shs,
00202 dir,
00203 navdir,
00204 iSetup);
00205 for (std::vector<TrajectorySeed*>::const_iterator iSeed = seeds.begin(); iSeed != seeds.end(); iSeed++){
00206 if (!*iSeed) {edm::LogError("CtfSpecialSeedGenerator")<<"a seed pointer is null. skipping.";continue;}
00207 if (postCheck(**iSeed)){
00208 output.push_back(**iSeed);
00209 }
00210 delete *iSeed;
00211 edm::LogVerbatim("CtfSpecialSeedGenerator") << "Seed built";
00212 }
00213 }
00214 }
00215 if ((theMaxSeeds > 0) && (output.size() > size_t(theMaxSeeds))) {
00216 edm::LogWarning("TooManySeeds") << "Too many seeds ("<< output.size() <<"), bailing out.\n";
00217 output.clear();
00218 return false;
00219 }
00220 return true;
00221 }
00222
00223 bool CtfSpecialSeedGenerator::preliminaryCheck(const SeedingHitSet& shs, const edm::EventSetup &es ){
00224
00225 edm::ESHandle<TrackerTopology> tTopo;
00226 es.get<IdealGeometryRecord>().get(tTopo);
00227
00228 std::vector<std::pair<unsigned int, unsigned int> > vSubdetLayer;
00229
00230 bool checkHitsAtPositiveY = conf_.getParameter<bool>("SeedsFromPositiveY");
00231
00232 bool checkHitsAtNegativeY = conf_.getParameter<bool>("SeedsFromNegativeY");
00233
00234 bool checkHitsOnDifferentLayers = conf_.getParameter<bool>("CheckHitsAreOnDifferentLayers");
00235 unsigned int nHits = shs.size();
00236 for (unsigned int iHit=0; iHit < nHits; ++iHit) {
00237
00238 const TrackingRecHit * trh = shs[iHit]->hit();
00239 TransientTrackingRecHit::RecHitPointer recHit = theBuilder->build(trh);
00240 GlobalPoint hitPos = recHit->globalPosition();
00241
00242
00243 if (checkHitsAtPositiveY){ if (hitPos.y() < 0) return false;}
00244
00245 if (checkHitsAtNegativeY){ if (hitPos.y() > 0) return false;}
00246
00247
00248
00249 unsigned int subid=(*trh).geographicalId().subdetId();
00250 unsigned int layer = tTopo->layer( (*trh).geographicalId());
00251 std::vector<std::pair<unsigned int, unsigned int> >::const_iterator iter;
00252
00253 if (checkHitsOnDifferentLayers){
00254
00255 for (iter = vSubdetLayer.begin(); iter != vSubdetLayer.end(); iter++){
00256 if (iter->first == subid && iter->second == layer) return false;
00257 }
00258
00259
00260
00261
00262
00263 }
00264
00265 vSubdetLayer.push_back(std::make_pair(subid, layer));
00266 }
00267 return true;
00268 }
00269
00270
00271 bool CtfSpecialSeedGenerator::postCheck(const TrajectorySeed& seed){
00272 if (!useScintillatorsConstraint) return true;
00273
00274 PTrajectoryStateOnDet pstate = seed.startingState();
00275 TrajectoryStateOnSurface theTSOS = trajectoryStateTransform::transientState(pstate,
00276 &(theTracker->idToDet(DetId(pstate.detId()))->surface()),
00277 &(*theMagfield));
00278 FreeTrajectoryState* state = theTSOS.freeState();
00279 StraightLinePlaneCrossing planeCrossingLower( Basic3DVector<float>(state->position()),
00280 Basic3DVector<float>(state->momentum()),
00281 alongMomentum);
00282 StraightLinePlaneCrossing planeCrossingUpper( Basic3DVector<float>(state->position()),
00283 Basic3DVector<float>(state->momentum()),
00284 oppositeToMomentum);
00285 std::pair<bool,StraightLinePlaneCrossing::PositionType> positionUpper =
00286 planeCrossingUpper.position(*upperScintillator);
00287 std::pair<bool,StraightLinePlaneCrossing::PositionType> positionLower =
00288 planeCrossingLower.position(*lowerScintillator);
00289 if (!(positionUpper.first && positionLower.first)) {
00290 edm::LogVerbatim("CtfSpecialSeedGenerator::checkDirection")
00291 << "Scintillator plane not crossed";
00292 return false;
00293 }
00294 LocalPoint positionUpperLocal = upperScintillator->toLocal((GlobalPoint)(positionUpper.second));
00295 LocalPoint positionLowerLocal = lowerScintillator->toLocal((GlobalPoint)(positionLower.second));
00296 if (upperScintillator->bounds().inside(positionUpperLocal) &&
00297 lowerScintillator->bounds().inside(positionLowerLocal)) {
00298 edm::LogVerbatim("CtfSpecialSeedGenerator::checkDirection")
00299 << "position on Upper scintillator "
00300 << positionUpper.second;
00301 edm::LogVerbatim("CtfSpecialSeedGenerator::checkDirection")
00302 << "position on Lower scintillator "
00303 << positionLower.second;
00304
00305 return true;
00306 }
00307 edm::LogVerbatim("CtfSpecialSeedGenerator::checkDirection")
00308 << "scintillator not crossed in bounds: position on Upper scintillator "
00309 << positionUpper.second << " position on Lower scintillator " << positionLower.second;
00310 return false;
00311 }
00312
00313