CMS 3D CMS Logo

Public Member Functions | Private Member Functions | Private Attributes

SimpleCosmicBONSeeder Class Reference

#include <SimpleCosmicBONSeeder.h>

Inheritance diagram for SimpleCosmicBONSeeder:
edm::EDProducer edm::ProducerBase edm::ProductRegistryHelper

List of all members.

Public Member Functions

void done ()
bool goodTriplet (const GlobalPoint &inner, const GlobalPoint &middle, const GlobalPoint &outer, const double &minRho) const
void init (const edm::EventSetup &c)
std::pair< GlobalVector, int > pqFromHelixFit (const GlobalPoint &inner, const GlobalPoint &middle, const GlobalPoint &outer, const edm::EventSetup &iSetup) const
virtual void produce (edm::Event &e, const edm::EventSetup &c)
bool seeds (TrajectorySeedCollection &output, const edm::EventSetup &iSetup)
 SimpleCosmicBONSeeder (const edm::ParameterSet &conf)
bool triplets (const edm::Event &e, const edm::EventSetup &c)
virtual ~SimpleCosmicBONSeeder ()

Private Member Functions

bool checkCharge (const TrackingRecHit *hit) const
bool checkCharge (const SiStripRecHit2D &hit, int subdetid) const
void checkNoisyModules (const std::vector< TransientTrackingRecHit::RecHitPointer > &hits, std::vector< bool > &oks) const

Private Attributes

std::string builderName
std::vector< int32_t > chargeThresholds_
ClusterChecker check_
bool checkCharge_
bool checkMaxHitsPerModule_
edm::ParameterSet conf_
int goodHitsPerSeed_
uint32_t helixVerbosity_
OrderedHitTriplets hitTriplets
std::vector< std::string > layerTripletNames_
edm::ESHandle< MagneticFieldmagfield
bool matchedRecHitUsesAnd_
std::vector< int32_t > maxHitsPerModule_
int32_t maxSeeds_
int32_t maxTriplets_
bool negativeYOnly
double pMin_
bool positiveYOnly
GlobalTrackingRegion region_
double rescaleError_
bool seedOnMiddle_
uint32_t seedVerbosity_
SeedingLayerSetsBuilder theLsb
PropagatorWithMaterialthePropagatorAl
PropagatorWithMaterialthePropagatorOp
KFUpdatortheUpdator
edm::ESHandle< TrackerGeometrytracker
uint32_t tripletsVerbosity_
edm::ESHandle
< TransientTrackingRecHitBuilder
TTTRHBuilder
bool writeTriplets_

Detailed Description

Definition at line 38 of file SimpleCosmicBONSeeder.h.


Constructor & Destructor Documentation

SimpleCosmicBONSeeder::SimpleCosmicBONSeeder ( const edm::ParameterSet conf) [explicit]

Definition at line 16 of file SimpleCosmicBONSeeder.cc.

References builderName, chargeThresholds_, checkCharge_, checkMaxHitsPerModule_, conf_, edm::ParameterSet::existsAs(), edm::ParameterSet::getParameter(), goodHitsPerSeed_, layerTripletNames_, matchedRecHitUsesAnd_, max(), maxHitsPerModule_, negativeYOnly, pMin_, positiveYOnly, ptmin, region_, sistripsummary::TEC, sistripsummary::TIB, sistripsummary::TID, sistripsummary::TOB, and writeTriplets_.

                                                                        : 
  conf_(conf),
  theLsb(conf.getParameter<edm::ParameterSet>("TripletsPSet")),
  writeTriplets_(conf.getParameter<bool>("writeTriplets")),
  seedOnMiddle_(conf.existsAs<bool>("seedOnMiddle") ? conf.getParameter<bool>("seedOnMiddle") : false),
  rescaleError_(conf.existsAs<double>("rescaleError") ? conf.getParameter<double>("rescaleError") : 1.0),
  tripletsVerbosity_(conf.getParameter<edm::ParameterSet>("TripletsPSet").getUntrackedParameter<uint32_t>("debugLevel",0)),
  seedVerbosity_(conf.getUntrackedParameter<uint32_t>("seedDebugLevel",0)),
  helixVerbosity_(conf.getUntrackedParameter<uint32_t>("helixDebugLevel",0)),
  check_(conf.getParameter<edm::ParameterSet>("ClusterCheckPSet")),
  maxTriplets_(conf.getParameter<int32_t>("maxTriplets")),
  maxSeeds_(conf.getParameter<int32_t>("maxSeeds"))
{
  edm::ParameterSet regionConf = conf_.getParameter<edm::ParameterSet>("RegionPSet");
  float ptmin        = regionConf.getParameter<double>("ptMin");
  float originradius = regionConf.getParameter<double>("originRadius");
  float halflength   = regionConf.getParameter<double>("originHalfLength");
  float originz      = regionConf.getParameter<double>("originZPosition");
  region_ = GlobalTrackingRegion(ptmin, originradius, halflength, originz);
  pMin_   = regionConf.getParameter<double>("pMin");

  builderName = conf_.getParameter<std::string>("TTRHBuilder");   

  //***top-bottom
  positiveYOnly=conf_.getParameter<bool>("PositiveYOnly");
  negativeYOnly=conf_.getParameter<bool>("NegativeYOnly");
  //***

  produces<TrajectorySeedCollection>();
  if (writeTriplets_) produces<edm::OwnVector<TrackingRecHit> >("cosmicTriplets");

  layerTripletNames_ = conf_.getParameter<edm::ParameterSet>("TripletsPSet").getParameter<std::vector<std::string> >("layerList");

  if (conf.existsAs<edm::ParameterSet>("ClusterChargeCheck")) {
      edm::ParameterSet cccc = conf.getParameter<edm::ParameterSet>("ClusterChargeCheck");
      checkCharge_          = cccc.getParameter<bool>("checkCharge");
      matchedRecHitUsesAnd_ = cccc.getParameter<bool>("matchedRecHitsUseAnd");
      chargeThresholds_.resize(7,0);
      edm::ParameterSet ccct = cccc.getParameter<edm::ParameterSet>("Thresholds");
      chargeThresholds_[StripSubdetector::TIB] = ccct.getParameter<int32_t>("TIB");
      chargeThresholds_[StripSubdetector::TID] = ccct.getParameter<int32_t>("TID");
      chargeThresholds_[StripSubdetector::TOB] = ccct.getParameter<int32_t>("TOB");
      chargeThresholds_[StripSubdetector::TEC] = ccct.getParameter<int32_t>("TEC");
  } else {
      checkCharge_ = false;
  }
  if (conf.existsAs<edm::ParameterSet>("HitsPerModuleCheck")) {
      edm::ParameterSet hpmcc = conf.getParameter<edm::ParameterSet>("HitsPerModuleCheck");
      checkMaxHitsPerModule_  = hpmcc.getParameter<bool>("checkHitsPerModule");
      maxHitsPerModule_.resize(7,std::numeric_limits<int32_t>::max());
      edm::ParameterSet hpmct = hpmcc.getParameter<edm::ParameterSet>("Thresholds");
      maxHitsPerModule_[StripSubdetector::TIB] = hpmct.getParameter<int32_t>("TIB");
      maxHitsPerModule_[StripSubdetector::TID] = hpmct.getParameter<int32_t>("TID");
      maxHitsPerModule_[StripSubdetector::TOB] = hpmct.getParameter<int32_t>("TOB");
      maxHitsPerModule_[StripSubdetector::TEC] = hpmct.getParameter<int32_t>("TEC");
  } else {
      checkMaxHitsPerModule_ = false;
  }
  if (checkCharge_ || checkMaxHitsPerModule_) {
      goodHitsPerSeed_ = conf.getParameter<int32_t>("minimumGoodHitsInSeed"); 
  } else {
      goodHitsPerSeed_ = 0;
  }

}
virtual SimpleCosmicBONSeeder::~SimpleCosmicBONSeeder ( ) [inline, virtual]

Definition at line 44 of file SimpleCosmicBONSeeder.h.

{}

Member Function Documentation

bool SimpleCosmicBONSeeder::checkCharge ( const TrackingRecHit hit) const [private]

Definition at line 305 of file SimpleCosmicBONSeeder.cc.

References cond::rpcobgas::detid, TrackingRecHit::geographicalId(), matchedRecHitUsesAnd_, and align::Tracker.

Referenced by triplets().

                                                                       {
    DetId detid(hit->geographicalId());
    if (detid.det() != DetId::Tracker) return false; // should not happen
    int subdet = detid.subdetId(); 
    if (subdet < 3) { // pixel
        return true;
    } else {
        if (typeid(*hit) == typeid(SiStripMatchedRecHit2D)) {
            const SiStripMatchedRecHit2D *mhit = static_cast<const SiStripMatchedRecHit2D *>(hit);
            if (matchedRecHitUsesAnd_) {
                return checkCharge(mhit->monoHit(), subdet) && checkCharge(mhit->stereoHit(), subdet);
            } else {
                return checkCharge(mhit->monoHit(), subdet) || checkCharge(mhit->stereoHit(), subdet);
            }
        } else if (typeid(*hit) == typeid(SiStripRecHit2D)) {
            return checkCharge(static_cast<const SiStripRecHit2D &>(*hit), subdet);
        } else {
            return true;
        }
    }
}
bool SimpleCosmicBONSeeder::checkCharge ( const SiStripRecHit2D &  hit,
int  subdetid 
) const [private]

Definition at line 328 of file SimpleCosmicBONSeeder.cc.

References SiStripCluster::amplitudes(), dtNoiseDBValidation_cfg::cerr, DeDxDiscriminatorTools::charge(), chargeThresholds_, SiStripCluster::firstStrip(), and tripletsVerbosity_.

                                                                                      {
    const SiStripCluster *clust = (hit.cluster().isNonnull() ?  hit.cluster().get() : hit.cluster_regional().get());
    int charge = std::accumulate(clust->amplitudes().begin(), clust->amplitudes().end(), int(0));
    if (tripletsVerbosity_ > 1) {
        std::cerr << "Hit on " << subdetid << ", charge = " << charge << ", threshold = " << chargeThresholds_[subdetid] 
                  << ", detid = " <<  hit.geographicalId().rawId() << ", firstStrip = " << clust->firstStrip() << std::endl;
    } else if ((tripletsVerbosity_ == 1) && (charge < chargeThresholds_[subdetid])) {
        std::cerr << "Hit on " << subdetid << ", charge = " << charge << ", threshold = " << chargeThresholds_[subdetid] 
                  << ", detid = " <<  hit.geographicalId().rawId() << ", firstStrip = " << clust->firstStrip() << std::endl;
    }
    return charge > chargeThresholds_[subdetid];
}
void SimpleCosmicBONSeeder::checkNoisyModules ( const std::vector< TransientTrackingRecHit::RecHitPointer > &  hits,
std::vector< bool > &  oks 
) const [private]

Definition at line 341 of file SimpleCosmicBONSeeder.cc.

References dtNoiseDBValidation_cfg::cerr, end, lumiContext::fill, maxHitsPerModule_, min, convertSQLiteXML::ok, DetId::rawId(), errorMatrix2Lands_multiChannel::start, DetId::subdetId(), and tripletsVerbosity_.

Referenced by triplets().

                                                                                                                                     {
    typedef TransientTrackingRecHit::RecHitPointer TTRH;
    std::vector<TTRH>::const_iterator it = hits.begin(),  start = it,   end = hits.end();
    std::vector<bool>::iterator       ok = oks.begin(), okStart = ok;
    while (start < end) {
        DetId lastid = (*start)->geographicalId();
        for (it = start + 1; (it < end) && ((*it)->geographicalId() == lastid); ++it) {
            ++ok;
        }
        if ( (it - start) > maxHitsPerModule_[lastid.subdetId()] ) { 
            if (tripletsVerbosity_ > 0) {
                std::cerr << "SimpleCosmicBONSeeder: Marking noisy module " << lastid.rawId() << ", it has " << (it-start) << " rechits"
                          << " (threshold is " << maxHitsPerModule_[lastid.subdetId()] << ")" << std::endl;
            }
            std::fill(okStart,ok,false);
        } else if (tripletsVerbosity_ > 0) {
            if ( (it - start) > std::min(4,maxHitsPerModule_[lastid.subdetId()]/4) ) {
                std::cerr << "SimpleCosmicBONSeeder: Not marking noisy module " << lastid.rawId() << ", it has " << (it-start) << " rechits"
                          << " (threshold is " << maxHitsPerModule_[lastid.subdetId()] << ")" << std::endl;
            }
        }
        start = it; okStart = ok;
    }
}
void SimpleCosmicBONSeeder::done ( )

Definition at line 604 of file SimpleCosmicBONSeeder.cc.

References thePropagatorAl, thePropagatorOp, and theUpdator.

Referenced by produce().

                                {
  delete thePropagatorAl;
  delete thePropagatorOp;
  delete theUpdator;
}
bool SimpleCosmicBONSeeder::goodTriplet ( const GlobalPoint inner,
const GlobalPoint middle,
const GlobalPoint outer,
const double &  minRho 
) const

Definition at line 366 of file SimpleCosmicBONSeeder.cc.

References gather_cfg::cout, FastCircle::rho(), tripletsVerbosity_, PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by triplets().

                                                                                                                                                   {
    float dyOM = outer.y() - middle.y(), dyIM = inner.y() - middle.y();
    if ((dyOM * dyIM > 0) && (fabs(dyOM)>10) && (fabs(dyIM)>10)) {
        if (tripletsVerbosity_ > 2) std::cout << "  fail for non coherent dy" << std::endl;
        return false;
    }
    float dzOM = outer.z() - middle.z(), dzIM = inner.z() - middle.z();
    if ((dzOM * dzIM > 0) && (fabs(dzOM)>50) && (fabs(dzIM)>50)) {
        if (tripletsVerbosity_ > 2) std::cout << "  fail for non coherent dz" << std::endl;
        return false;
    }
    if (minRho > 0) {
        FastCircle theCircle(inner,middle,outer);
        if (theCircle.rho() < minRho) {
            if (tripletsVerbosity_ > 2) std::cout << "  fail for pt cut" << std::endl;
            return false;
        }
    }
    return true;
}
void SimpleCosmicBONSeeder::init ( const edm::EventSetup c)
std::pair< GlobalVector, int > SimpleCosmicBONSeeder::pqFromHelixFit ( const GlobalPoint inner,
const GlobalPoint middle,
const GlobalPoint outer,
const edm::EventSetup iSetup 
) const

Definition at line 388 of file SimpleCosmicBONSeeder.cc.

References abs, GlobalTrajectoryParameters::charge(), gather_cfg::cout, helixVerbosity_, magfield, GlobalTrajectoryParameters::momentum(), FreeTrajectoryState::parameters(), FastCircle::rho(), rho, FastHelix::stateAtVertex(), PV3DBase< T, PVType, FrameType >::x(), FastCircle::x0(), PV3DBase< T, PVType, FrameType >::y(), FastCircle::y0(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by seeds(), and triplets().

                                                                                            {
    if (helixVerbosity_ > 0) {
        std::cout << "DEBUG PZ =====" << std::endl;
        FastHelix helix(inner,middle,outer,iSetup);
        GlobalVector gv=helix.stateAtVertex().parameters().momentum(); // status on inner hit
        std::cout << "FastHelix P = " << gv   << "\n";
        std::cout << "FastHelix Q = " << helix.stateAtVertex().parameters().charge() << "\n";
    }

    // My attempt (with different approx from FastHelix)
    // 1) fit the circle
    FastCircle theCircle(inner,middle,outer);
    double rho = theCircle.rho();

    // 2) Get the PT
    GlobalVector tesla = magfield->inTesla(middle);
    double pt = 0.01 * rho * (0.3*tesla.z());

    // 3) Get the PX,PY at OUTER hit (VERTEX)
    double dx1 = outer.x()-theCircle.x0();
    double dy1 = outer.y()-theCircle.y0();
    double py = pt*dx1/rho, px = -pt*dy1/rho;
    if(px*(middle.x() - outer.x()) + py*(middle.y() - outer.y()) < 0.) {
        px *= -1.; py *= -1.;
    }

    // 4) Get the PZ through pz = pT*(dz/d(R*phi)))
    double dz = inner.z() - outer.z();
    double sinphi = ( dx1*(inner.y()-theCircle.y0()) - dy1*(inner.x()-theCircle.x0())) / (rho * rho);
    double dphi = std::abs(std::asin(sinphi));
    double pz = pt * dz / (dphi * rho); 

    int myq = ((theCircle.x0()*py - theCircle.y0()*px) / tesla.z()) > 0. ?  +1 : -1;
    
    std::pair<GlobalVector,int> mypq(GlobalVector(px,py,pz),myq);

    if (helixVerbosity_ > 1) {
        std::cout << "Gio: pt = " << pt << std::endl;
        std::cout << "Gio: dz = " << dz << ", sinphi = " << sinphi << ", dphi = " << dphi << ", dz/drphi = " << (dz/dphi/rho) << std::endl;
    }
    if (helixVerbosity_ > 0) {
        std::cout << "Gio's fit P = " << mypq.first << "\n";
        std::cout << "Gio's fit Q = " << myq  << "\n";
    }

    return mypq;
}
void SimpleCosmicBONSeeder::produce ( edm::Event e,
const edm::EventSetup c 
) [virtual]

Implements edm::EDProducer.

Definition at line 83 of file SimpleCosmicBONSeeder.cc.

References check_, TrackingRecHit::clone(), done(), edm::EventSetup::get(), hitTriplets, init(), magfield, convertSQLitetoXML_cfg::output, edm::Event::put(), seeds(), ClusterChecker::tooManyClusters(), triplets(), and writeTriplets_.

{
  std::auto_ptr<TrajectorySeedCollection> output(new TrajectorySeedCollection());
  std::auto_ptr<edm::OwnVector<TrackingRecHit> > outtriplets(new edm::OwnVector<TrackingRecHit>());

  es.get<IdealMagneticFieldRecord>().get(magfield);
  if (magfield->inTesla(GlobalPoint(0,0,0)).mag() > 0.01) {
     size_t clustsOrZero = check_.tooManyClusters(ev);
     if (clustsOrZero) {
         edm::LogError("TooManyClusters") << "Found too many clusters (" << clustsOrZero << "), bailing out.\n";
     } else {
         init(es);
         bool tripletsOk = triplets(ev,es);
         if (tripletsOk) {

             bool seedsOk    = seeds(*output,es);
             if (!seedsOk) { }

             if (writeTriplets_) {
                 for (OrderedHitTriplets::const_iterator it = hitTriplets.begin(); it != hitTriplets.end(); ++it) {
                     const TrackingRecHit * hit1 = it->inner()->hit();    
                     const TrackingRecHit * hit2 = it->middle()->hit();
                     const TrackingRecHit * hit3 = it->outer()->hit();
                     outtriplets->push_back(hit1->clone());
                     outtriplets->push_back(hit2->clone());
                     outtriplets->push_back(hit3->clone());
                 }
             }
         }
         done();
     }
  }

  if (writeTriplets_) {
      ev.put(outtriplets, "cosmicTriplets");
  }
  ev.put(output);
}
bool SimpleCosmicBONSeeder::seeds ( TrajectorySeedCollection output,
const edm::EventSetup iSetup 
)

Definition at line 437 of file SimpleCosmicBONSeeder.cc.

References alongMomentum, FreeTrajectoryState::cartesianError(), gather_cfg::cout, cmsPerfPublish::fail(), hitTriplets, SurfaceOrientation::inner, OrderedHitTriplet::inner(), edm::detail::isnan(), magfield, CartesianTrajectoryError::matrix(), maxSeeds_, OrderedHitTriplet::middle(), oppositeToMomentum, SurfaceOrientation::outer, OrderedHitTriplet::outer(), PV3DBase< T, PVType, FrameType >::perp(), trajectoryStateTransform::persistentState(), pMin_, pqFromHelixFit(), Propagator::propagate(), LargeD0_PixelPairStep_cff::propagator, TrackingRegionBase::ptMin(), edm::OwnVector< T, P >::push_back(), region_, FreeTrajectoryState::rescaleError(), rescaleError_, seedOnMiddle_, seedVerbosity_, OrderedHitTriplets::size(), mathSSE::sqrt(), swap(), thePropagatorAl, thePropagatorOp, theUpdator, tracker, KFUpdator::update(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by produce().

{
    typedef TrajectoryStateOnSurface TSOS;
    
    for (size_t it=0;it<hitTriplets.size();it++){
        const OrderedHitTriplet &trip = hitTriplets[it];

        GlobalPoint inner = tracker->idToDet((*(trip.inner())).geographicalId())->surface().
            toGlobal((*(trip.inner())).localPosition());

        GlobalPoint middle = tracker->idToDet((*(trip.middle())).geographicalId())->surface().
            toGlobal((*(trip.middle())).localPosition());

        GlobalPoint outer = tracker->idToDet((*(trip.outer())).geographicalId())->surface().
            toGlobal((*(trip.outer())).localPosition());   

        if (seedVerbosity_ > 1)
            std::cout << "Processing triplet " << it << ": " << inner << " + " << middle << " + " << outer << std::endl;

        if ( (outer.y()-inner.y())*outer.y() < 0 ) {
            std::swap(inner,outer);
            std::swap(const_cast<TransientTrackingRecHit::ConstRecHitPointer &>(trip.inner()),
                      const_cast<TransientTrackingRecHit::ConstRecHitPointer &>(trip.outer()) );

//            std::swap(const_cast<ctfseeding::SeedingHit &>(trip.inner()), 
//                      const_cast<ctfseeding::SeedingHit &>(trip.outer()) );
            if (seedVerbosity_ > 1) {
                std::cout << "The seed was going away from CMS! swapped in <-> out" << std::endl;
                std::cout << "Processing swapped triplet " << it << ": " << inner << " + " << middle << " + " << outer << std::endl;
            }
        }

        // First use FastHelix out of the box
        std::pair<GlobalVector,int> pq = pqFromHelixFit(inner,middle,outer,iSetup);
        GlobalVector gv = pq.first;
        float        ch = pq.second; 
        float Mom = sqrt( gv.x()*gv.x() + gv.y()*gv.y() + gv.z()*gv.z() ); 

        if(Mom > 10000 || isnan(Mom))  { 
            if (seedVerbosity_ > 1)
                std::cout << "Processing triplet " << it << ": fail for momentum." << std::endl; 
            continue;
        }

        if (gv.perp() < region_.ptMin()) {
            if (seedVerbosity_ > 1)
                std::cout << "Processing triplet " << it << ": fail for pt = " << gv.perp() << " < ptMin = " << region_.ptMin() << std::endl; 
            continue;
        }

        const Propagator * propagator = 0;  
        if((outer.y()-inner.y())>0){
            if (seedVerbosity_ > 1)
                std::cout << "Processing triplet " << it << ":  downgoing." << std::endl; 
            propagator = thePropagatorAl;
        } else {
            gv = -1*gv; ch = -1.*ch;                        
            propagator = thePropagatorOp;
            if (seedVerbosity_ > 1)
                std::cout << "Processing triplet " << it << ":  upgoing." << std::endl; 
        }

        if (seedVerbosity_ > 1) {
            if (( gv.z() * (outer.z()-inner.z()) > 0 ) && ( fabs(outer.z()-inner.z()) > 5) && (fabs(gv.z()) > .01))  {
                std::cout << "ORRORE: outer.z()-inner.z() = " << (outer.z()-inner.z()) << ", gv.z() = " << gv.z() << std::endl;
            }
        }

        GlobalTrajectoryParameters Gtp(outer,
                gv,int(ch), 
                &(*magfield));
        FreeTrajectoryState CosmicSeed(Gtp,
                CurvilinearTrajectoryError(AlgebraicSymMatrix55(AlgebraicMatrixID())));  
        CosmicSeed.rescaleError(100);
        if (seedVerbosity_ > 2) {
            std::cout << "Processing triplet " << it << ". start from " << std::endl;
            std::cout << "    X  = " << outer << ", P = " << gv << std::endl;
            std::cout << "    Cartesian error (X,P) = \n" << CosmicSeed.cartesianError().matrix() << std::endl;
        }
       
        edm::OwnVector<TrackingRecHit> hits;
        OrderedHitTriplet seedHits(trip.outer(),trip.middle(),trip.inner());
        TSOS propagated, updated;
        bool fail = false;
        for (size_t ih = 0; ih < 3; ++ih) {
            if ((ih == 2) && seedOnMiddle_) {
                if (seedVerbosity_ > 2) 
                    std::cout << "Stopping at middle hit, as requested." << std::endl;
                break;
            }
            if (seedVerbosity_ > 2)
                std::cout << "Processing triplet " << it << ", hit " << ih << "." << std::endl;
            if (ih == 0) {
                propagated = propagator->propagate(CosmicSeed, tracker->idToDet((*seedHits[ih]).geographicalId())->surface());
            } else {
                propagated = propagator->propagate(updated, tracker->idToDet((*seedHits[ih]).geographicalId())->surface());
            }
            if (!propagated.isValid()) {
                if (seedVerbosity_ > 1)
                    std::cout << "Processing triplet " << it << ", hit " << ih << ": failed propagation." << std::endl;
                fail = true; break;
            } else {
                if (seedVerbosity_ > 2)
                    std::cout << "Processing triplet " << it << ", hit " << ih << ": propagated state = " << propagated;
            }
            const TransientTrackingRecHit::ConstRecHitPointer & tthp   = seedHits[ih];
            TransientTrackingRecHit::RecHitPointer              newtth = tthp->clone(propagated);
            hits.push_back(newtth->hit()->clone());
            updated = theUpdator->update(propagated, *newtth);
            if (!updated.isValid()) {
                if (seedVerbosity_ > 1)
                    std::cout << "Processing triplet " << it << ", hit " << ih << ": failed update." << std::endl;
                fail = true; break;
            } else {
                if (seedVerbosity_ > 2)
                    std::cout << "Processing triplet " << it << ", hit " << ih << ": updated state = " << updated;
            }
        }
        if (!fail && updated.isValid() && (updated.globalMomentum().perp() < region_.ptMin())) {
            if (seedVerbosity_ > 1)
                std::cout << "Processing triplet " << it << 
                             ": failed for final pt " << updated.globalMomentum().perp() << " < " << region_.ptMin() << std::endl;
            fail = true;
        }
        if (!fail && updated.isValid() && (updated.globalMomentum().mag() < pMin_)) {
            if (seedVerbosity_ > 1)
                std::cout << "Processing triplet " << it << 
                             ": failed for final p " << updated.globalMomentum().perp() << " < " << pMin_ << std::endl;
            fail = true;
        }
        if (!fail) {
            if (rescaleError_ != 1.0) {
                if (seedVerbosity_ > 2) {
                    std::cout << "Processing triplet " << it << ", rescale error by " << rescaleError_ << ": state BEFORE rescaling " << updated;
                    std::cout << "    Cartesian error (X,P) before rescaling= \n" << updated.cartesianError().matrix() << std::endl;
                }
                updated.rescaleError(rescaleError_);
            }
            if (seedVerbosity_ > 0) {
                std::cout << "Processed  triplet " << it << ": success (saved as #"<<output.size()<<") : " 
                        << inner << " + " << middle << " + " << outer << std::endl;
                std::cout << "    pt = " << updated.globalMomentum().perp() <<
                             "    eta = " << updated.globalMomentum().eta() << 
                             "    phi = " << updated.globalMomentum().phi() <<
                             "    ch = " << updated.charge() << std::endl;
                if (seedVerbosity_ > 1) {
                    std::cout << "    State:" << updated;
                } else {
                    std::cout << "    X  = " << updated.globalPosition() << ", P = " << updated.globalMomentum() << std::endl;
                }
                std::cout << "    Cartesian error (X,P) = \n" << updated.cartesianError().matrix() << std::endl;
            }
            
            PTrajectoryStateOnDet const &  PTraj = trajectoryStateTransform::persistentState(updated, 
                                                            (*(seedOnMiddle_ ? trip.middle() : trip.inner())).geographicalId().rawId());
            output.push_back(TrajectorySeed(PTraj,hits,
                                                ( (outer.y()-inner.y()>0) ? alongMomentum : oppositeToMomentum) ));
            if ((maxSeeds_ > 0) && (output.size() > size_t(maxSeeds_))) { 
                output.clear(); 
                edm::LogError("TooManySeeds") << "Found too many seeds, bailing out.\n";
                return false;
            }
        }
    }
    return true;
}
bool SimpleCosmicBONSeeder::triplets ( const edm::Event e,
const edm::EventSetup c 
)

ctfseeding SeedinHits and their iterators

Transient Tracking RecHits

Checks on the cluster charge and on noisy modules

Now actually filling in the charges for all the clusters

Definition at line 168 of file SimpleCosmicBONSeeder.cc.

References newFWLiteAna::build, checkCharge(), checkCharge_, checkMaxHitsPerModule_, checkNoisyModules(), gather_cfg::cout, Exception, goodHitsPerSeed_, goodTriplet(), helixVerbosity_, GlobalTrackingRegion::hits(), hitTriplets, UserOptions_cff::idx, SeedingLayerSetsBuilder::layers(), layerTripletNames_, python::rootplot::utilities::ls(), magfield, maxTriplets_, negativeYOnly, positiveYOnly, pqFromHelixFit(), TrackingRegionBase::ptMin(), region_, OrderedHitTriplets::size(), python::multivaluedict::sort(), theLsb, tripletsVerbosity_, and PV3DBase< T, PVType, FrameType >::y().

Referenced by produce().

                                                                               {
    using namespace ctfseeding;

    hitTriplets.clear();
    hitTriplets.reserve(0);
    SeedingLayerSets lss = theLsb.layers(es);
    SeedingLayerSets::const_iterator iLss;

    double minRho = region_.ptMin() / ( 0.003 * magfield->inTesla(GlobalPoint(0,0,0)).z() );

    for (iLss = lss.begin(); iLss != lss.end(); iLss++){
        SeedingLayers ls = *iLss;
        if (ls.size() != 3){
            throw cms::Exception("CtfSpecialSeedGenerator") << "You are using " << ls.size() <<" layers in set instead of 3 ";
        }

        std::vector<SeedingHit> innerHits  = region_.hits(e, es, &ls[0]);
        std::vector<SeedingHit> middleHits = region_.hits(e, es, &ls[1]);
        std::vector<SeedingHit> outerHits  = region_.hits(e, es, &ls[2]);
        std::vector<SeedingHit>::const_iterator iOuterHit,iMiddleHit,iInnerHit;

        if (tripletsVerbosity_ > 0) {
            std::cout << "GenericTripletGenerator iLss = " << layerTripletNames_[iLss - lss.begin()]
                    << " (" << (iLss - lss.begin()) << "): # = " 
                    << innerHits.size() << "/" << middleHits.size() << "/" << outerHits.size() << std::endl;
        }

        typedef TransientTrackingRecHit::RecHitPointer TTRH;
        std::vector<TTRH> innerTTRHs, middleTTRHs, outerTTRHs;

        std::vector<bool> innerOk( innerHits.size(),  true);
        std::vector<bool> middleOk(middleHits.size(), true);
        std::vector<bool> outerOk( outerHits.size(),  true);

        size_t sizBefore = hitTriplets.size();
        int idx = 0;
        for (iOuterHit = outerHits.begin(), idx = 0; iOuterHit != outerHits.end(); ++idx, ++iOuterHit){
            outerTTRHs.push_back(ls[2].hitBuilder()->build((**iOuterHit).hit()));
            if (checkCharge_ && !checkCharge(outerTTRHs.back()->hit())) outerOk[idx] = false;
        }
        for (iMiddleHit = middleHits.begin(), idx = 0; iMiddleHit != middleHits.end(); ++idx, ++iMiddleHit){
            middleTTRHs.push_back(ls[1].hitBuilder()->build((**iMiddleHit).hit()));
            if (checkCharge_ && !checkCharge(middleTTRHs.back()->hit())) middleOk[idx] = false;
        }
        for (iInnerHit = innerHits.begin(), idx = 0; iInnerHit != innerHits.end(); ++idx, ++iInnerHit){
            innerTTRHs.push_back(ls[0].hitBuilder()->build((**iInnerHit).hit()));
            if (checkCharge_ && !checkCharge(innerTTRHs.back()->hit())) innerOk[idx] = false;
        }
        if (checkMaxHitsPerModule_) {
            checkNoisyModules(innerTTRHs,  innerOk);
            checkNoisyModules(middleTTRHs, middleOk);
            checkNoisyModules(outerTTRHs,  outerOk);
        }

        for (iOuterHit = outerHits.begin(); iOuterHit != outerHits.end(); iOuterHit++){
            idx = iOuterHit - outerHits.begin();
            TTRH &       outerTTRH = outerTTRHs[idx];
            GlobalPoint  outerpos  = outerTTRH->globalPosition(); // this caches by itself
            bool         outerok   = outerOk[idx];
            if (outerok < goodHitsPerSeed_ - 2) {
                if (tripletsVerbosity_ > 2) 
                    std::cout << "Skipping at first hit: " << (outerok) << " < " << (goodHitsPerSeed_ - 2) << std::endl;
                continue; 
            }

            for (iMiddleHit = middleHits.begin(); iMiddleHit != middleHits.end(); iMiddleHit++){
                idx = iMiddleHit - middleHits.begin();
                TTRH &       middleTTRH = middleTTRHs[idx];
                GlobalPoint  middlepos  = middleTTRH->globalPosition(); // this caches by itself
                bool         middleok   = middleOk[idx];
                if (outerok+middleok < goodHitsPerSeed_ - 1) {
                    if (tripletsVerbosity_ > 2) 
                        std::cout << "Skipping at second hit: " << (outerok+middleok) << " < " << (goodHitsPerSeed_ - 1) << std::endl;
                    continue; 
                }

                for (iInnerHit = innerHits.begin(); iInnerHit != innerHits.end(); iInnerHit++){
                    idx = iInnerHit - innerHits.begin();
                    TTRH &       innerTTRH = innerTTRHs[idx];
                    GlobalPoint  innerpos  = innerTTRH->globalPosition(); // this caches by itself
                    bool         innerok   = innerOk[idx];
                    if (outerok+middleok+innerok < goodHitsPerSeed_) {
                        if (tripletsVerbosity_ > 2) 
                            std::cout << "Skipping at third hit: " << (outerok+middleok+innerok) << " < " << (goodHitsPerSeed_) << std::endl;
                        continue;
                    } 

                    //***top-bottom
                    if (positiveYOnly && (innerpos.y()<0 || middlepos.y()<0 || outerpos.y()<0
                                          || outerpos.y() < innerpos.y()
                                          ) ) continue;
                    if (negativeYOnly && (innerpos.y()>0 || middlepos.y()>0 || outerpos.y()>0
                                          || outerpos.y() > innerpos.y()
                                          ) ) continue;
                    //***
                    
                    if (tripletsVerbosity_ > 2) std::cout << "Trying seed with: " << innerpos << " + " << middlepos << " + " << outerpos << std::endl;
                    if (goodTriplet(innerpos,middlepos,outerpos,minRho)) {
                        OrderedHitTriplet oht(*iInnerHit,*iMiddleHit,*iOuterHit);
                        hitTriplets.push_back(oht);
                        if ((maxTriplets_ > 0) && (hitTriplets.size() > size_t(maxTriplets_))) {
                            hitTriplets.clear();                      // clear
                            //OrderedHitTriplets().swap(hitTriplets); // really clear   
                            edm::LogError("TooManyTriplets") << "Found too many triplets, bailing out.\n";
                            return false;
                        }
                        if (tripletsVerbosity_ > 3) {
                            std::cout << " accepted seed #" << (hitTriplets.size()-1) << " w/: " 
                                << innerpos << " + " << middlepos << " + " << outerpos << std::endl;
                        }
                        if (tripletsVerbosity_ == 2) {
                                std::cout << " good seed #" << (hitTriplets.size()-1) << " w/: "     
                                    << innerpos << " + " << middlepos << " + " << outerpos << std::endl;
                        }
                        if (tripletsVerbosity_ > 3 && (helixVerbosity_ > 0)) { // debug the momentum here too
                            pqFromHelixFit(innerpos,middlepos,outerpos,es); 
                        }
                    }
                }
            }
        }
        if ((tripletsVerbosity_ > 0) && (hitTriplets.size() > sizBefore)) {
            std::cout << "                        iLss = " << layerTripletNames_[iLss - lss.begin()]
                << " (" << (iLss - lss.begin()) << "): # = " 
                << innerHits.size() << "/" << middleHits.size() << "/" << outerHits.size() 
                << ": Found " << (hitTriplets.size() - sizBefore) << " seeds [running total: " << hitTriplets.size() << "]"
                << std::endl ;
        }

    }
    std::sort(hitTriplets.begin(),hitTriplets.end(),HigherInnerHit());
    return true;
}

Member Data Documentation

std::string SimpleCosmicBONSeeder::builderName [private]

Definition at line 60 of file SimpleCosmicBONSeeder.h.

Referenced by init(), and SimpleCosmicBONSeeder().

std::vector<int32_t> SimpleCosmicBONSeeder::chargeThresholds_ [private]

Definition at line 89 of file SimpleCosmicBONSeeder.h.

Referenced by checkCharge(), and SimpleCosmicBONSeeder().

Definition at line 81 of file SimpleCosmicBONSeeder.h.

Referenced by produce().

Definition at line 87 of file SimpleCosmicBONSeeder.h.

Referenced by SimpleCosmicBONSeeder(), and triplets().

Definition at line 90 of file SimpleCosmicBONSeeder.h.

Referenced by SimpleCosmicBONSeeder(), and triplets().

Definition at line 59 of file SimpleCosmicBONSeeder.h.

Referenced by SimpleCosmicBONSeeder().

Definition at line 86 of file SimpleCosmicBONSeeder.h.

Referenced by SimpleCosmicBONSeeder(), and triplets().

Definition at line 70 of file SimpleCosmicBONSeeder.h.

Referenced by pqFromHelixFit(), and triplets().

Definition at line 84 of file SimpleCosmicBONSeeder.h.

Referenced by produce(), seeds(), and triplets().

std::vector<std::string> SimpleCosmicBONSeeder::layerTripletNames_ [private]

Definition at line 71 of file SimpleCosmicBONSeeder.h.

Referenced by SimpleCosmicBONSeeder(), and triplets().

Definition at line 73 of file SimpleCosmicBONSeeder.h.

Referenced by init(), pqFromHelixFit(), produce(), seeds(), and triplets().

Definition at line 88 of file SimpleCosmicBONSeeder.h.

Referenced by checkCharge(), and SimpleCosmicBONSeeder().

std::vector<int32_t> SimpleCosmicBONSeeder::maxHitsPerModule_ [private]

Definition at line 91 of file SimpleCosmicBONSeeder.h.

Referenced by checkNoisyModules(), and SimpleCosmicBONSeeder().

Definition at line 82 of file SimpleCosmicBONSeeder.h.

Referenced by seeds().

Definition at line 82 of file SimpleCosmicBONSeeder.h.

Referenced by triplets().

Definition at line 98 of file SimpleCosmicBONSeeder.h.

Referenced by SimpleCosmicBONSeeder(), and triplets().

double SimpleCosmicBONSeeder::pMin_ [private]

Definition at line 64 of file SimpleCosmicBONSeeder.h.

Referenced by seeds(), and SimpleCosmicBONSeeder().

Definition at line 97 of file SimpleCosmicBONSeeder.h.

Referenced by SimpleCosmicBONSeeder(), and triplets().

Definition at line 63 of file SimpleCosmicBONSeeder.h.

Referenced by seeds(), SimpleCosmicBONSeeder(), and triplets().

Definition at line 68 of file SimpleCosmicBONSeeder.h.

Referenced by seeds().

Definition at line 67 of file SimpleCosmicBONSeeder.h.

Referenced by seeds().

Definition at line 70 of file SimpleCosmicBONSeeder.h.

Referenced by seeds().

Definition at line 62 of file SimpleCosmicBONSeeder.h.

Referenced by triplets().

Definition at line 77 of file SimpleCosmicBONSeeder.h.

Referenced by done(), init(), and seeds().

Definition at line 78 of file SimpleCosmicBONSeeder.h.

Referenced by done(), init(), and seeds().

Definition at line 76 of file SimpleCosmicBONSeeder.h.

Referenced by done(), init(), and seeds().

Definition at line 74 of file SimpleCosmicBONSeeder.h.

Referenced by init(), and seeds().

Definition at line 70 of file SimpleCosmicBONSeeder.h.

Referenced by checkCharge(), checkNoisyModules(), goodTriplet(), and triplets().

Definition at line 75 of file SimpleCosmicBONSeeder.h.

Referenced by init().

Definition at line 65 of file SimpleCosmicBONSeeder.h.

Referenced by produce(), and SimpleCosmicBONSeeder().