50 #include <boost/regex.hpp>
75 namespace reco {
namespace modules {
87 Rule(
const std::string &str) ;
165 static boost::regex rule(
"(keep|drop)\\s+([A-Z]+)(\\s+(\\d+))?");
171 if (!regex_match(str.c_str(),
match, rule)) {
172 throw cms::Exception(
"Configuration") <<
"Rule '" << str <<
"' not understood.\n";
175 std::cout<<
"*** Rule Command given to TrackerTrackHitFilter:\t"<<str<<std::endl;
180 keep_ = (strncmp(match[1].
first,
"keep",4 )== 0);
191 throw cms::Exception(
"Configuration") <<
"Detector '" << match[2].first <<
"' not understood. Should be PXB, PXE, TIB, TID, TOB, TEC.\n";
194 if (match[4].first != match[4].
second) {
195 layer_ = atoi(match[4].first);
218 static boost::regex rule(
"([A-Z]+)"
219 "\\s*(\\d+\\.*\\d*)?"
220 "\\s*(\\d+\\.*\\d*)?");
228 if (!regex_match(str.c_str(),
match, rule)) {
229 throw cms::Exception(
"Configuration") <<
"Rule for S to N cut '" << str <<
"' not understood.\n";
232 std::string match_0=(match[0].first,match[0].second);
233 match_1=(match[1].first,match[1].second);
234 match_2=(match[2].first,match[2].second);
235 match_3=(match[3].first,match[3].second);
241 for(cnt=0;cnt<6;cnt++ ){
242 subdet_ind[cnt]=-1.0;
247 std::string match_1a(match[1].
first,match[1].
second);
248 if (strncmp(match[1].first,
"ALL", 3) == 0) doALL=
true;
249 if (doALL || strncmp(match[1].first,
"PXB", 3) == 0) subdet_ind[0] = +1.0;
250 if (doALL || strncmp(match[1].first,
"PXE", 3) == 0) subdet_ind[1] = +1.0;
251 if (doALL || strncmp(match[1].first,
"TIB", 3) == 0) subdet_ind[2] = +1.0;
252 if (doALL || strncmp(match[1].first,
"TID", 3) == 0) subdet_ind[3] = +1.0;
253 if (doALL || strncmp(match[1].first,
"TOB", 3) == 0) subdet_ind[4] = +1.0;
254 if (doALL || strncmp(match[1].first,
"TEC", 3) == 0) subdet_ind[5] = +1.0;
256 for(cnt=0;cnt<6;cnt++ ){
257 if(subdet_ind[cnt]>0.0){
259 if (match[2].first != match[2].second) {
262 if (match[3].first != match[3].second ) {
269 bool correct_regex=
false;
270 for(cnt=0;cnt<6;cnt++ ){
274 if (!correct_regex) {
275 throw cms::Exception(
"Configuration") <<
"Detector '" << match_1a <<
"' not understood in parseStoN. Should be PXB, PXE, TIB, TID, TOB, TEC.\n";
283 src_(iConfig.getParameter<edm::InputTag>(
"src")),
284 minimumHits_(iConfig.getParameter<uint32_t>(
"minimumHits")),
291 detsToIgnore_( iConfig.getParameter<std::vector<uint32_t> >(
"detsToIgnore") ),
296 pxlTPLProbXY_(iConfig.getParameter<double>(
"PxlTemplateProbXYCut")),
297 pxlTPLProbXYQ_(iConfig.getParameter<double>(
"PxlTemplateProbXYChargeCut")),
298 pxlTPLqBin_(iConfig.getParameter<std::vector<int32_t> >(
"PxlTemplateqBinCut")),
300 tagOverlaps_( iConfig.getParameter<bool>(
"tagOverlaps") )
306 throw cms::Exception(
"Configuration") <<
"Inconsistent Configuration: you can't set both 'stripAllInvalidHits' and 'replaceWithInactiveHits' to true\n";
309 throw cms::Exception(
"Configuration") <<
"Wrong configuration of TrackerTrackHitFilter. You cannot apply the cut on the track angle without using Trajectories!\n";
312 throw cms::Exception(
"Configuration") <<
"Wrong configuration of TrackerTrackHitFilter. You cannot apply the cut on the corrected pixel cluster charge without using Trajectories!\n";
316 std::cout<<
"Warning from TrackerTrackHitFilter: vector with qBin cuts has size > 2. Additional items will be ignored."<<std::endl;
322 std::vector<std::string> str_rules = iConfig.
getParameter<std::vector<std::string> >(
"commands");
323 rules_.reserve(str_rules.size());
324 for (std::vector<std::string>::const_iterator it = str_rules.begin(), ed = str_rules.end(); it != ed; ++it) {
334 for(cnt=0;cnt<6;cnt++ ){
340 std::vector<std::string> str_StoNrules = iConfig.
getParameter<std::vector<std::string> >(
"StoNcommands");
341 for (std::vector<std::string>::const_iterator str_StoN = str_StoNrules.begin(); str_StoN != str_StoNrules.end(); ++str_StoN) {
345 std::cout<<
"Finished parsing S/N. Applying following cuts to subdets:";
346 for(cnt=0;cnt<6;cnt++ ){
360 produces<TrackCandidateCollection>();
387 else candcollsize=tracks->size();
390 output->reserve(candcollsize);
393 std::vector<TrackingRecHit *> hits;
402 const Track *trk = &(*tkref);
406 std::vector<TrackingRecHit *>::iterator
begin = hits.begin(),
end = hits.end();
410 while ( (begin !=
end) && ( (*begin)->isValid() ==
false ) ) ++
begin;
415 while ( (begin !=
end) && ( (*end)->isValid() ==
false ) ) --
end;
422 for( std::vector<TrackingRecHit *>::iterator ithit=begin;ithit!=
end;++ithit){
423 if( (*ithit)->isValid())nvalidhits++;
444 for (begin = hits.begin(),
end = hits.end(); begin !=
end; ++
begin) {
445 if (*begin)
delete *
begin;
455 for (std::vector<reco::Track>::const_iterator ittrk = tracks->begin(), edtrk = tracks->end(); ittrk != edtrk; ++ittrk) {
459 const Track *trk = &(*ittrk);
485 std::vector<TrackingRecHit *>::iterator
begin = hits.begin(),
end = hits.end();
489 while ( (begin !=
end) && ( (*begin)->isValid() ==
false ) ) ++
begin;
494 while ( (begin !=
end) && ( (*end)->isValid() ==
false ) ) --
end;
501 for( std::vector<TrackingRecHit *>::iterator ithit=begin;ithit!=
end;++ithit){
502 if( (*ithit)->isValid())nvalidhits++;
516 for (begin = hits.begin(),
end = hits.end(); begin !=
end; ++
begin) {
517 if (*begin)
delete *
begin;
534 if ( pdir ==
anyDirection )
throw cms::Exception(
"UnimplementedFeature") <<
"Cannot work with tracks that have 'anyDirecton' \n";
550 ownHits.
reserve(hitsEnd - hitsBegin);
551 for ( ; hitsBegin != hitsEnd; ++hitsBegin) {
572 if(hit->
isValid() && hit==0 && detid.
rawId()==0)
continue;
574 int verdict=
checkHit(iSetup,detid,hit);
577 hits.push_back(hit->
clone());
585 else if(verdict==-2) hits.push_back(hit->
clone());
586 else if(verdict==-1){
588 hits.push_back(hit->
clone());
601 std::vector<TrajectoryMeasurement> tmColl =itt->
measurements();
611 for(std::vector<TrajectoryMeasurement>::const_iterator itTrajMeas = tmColl.begin(); itTrajMeas!=tmColl.end(); itTrajMeas++){
615 if(hitpointer->isValid() && hitpointer->hit()==0){constrhits++;
continue;}
619 int verdict=
checkHit(iSetup,detid,hit);
646 if ( ( previousTM!=0 )&& (layer!=-1 )) {
648 for (std::vector<TrajectoryMeasurement>::const_iterator itmCompare =itTrajMeas-1;itmCompare >= tmColl.begin() && itmCompare > itTrajMeas - 4;--itmCompare){
650 DetId compareId = itmCompare->recHit()->geographicalId();
651 if ( subDet != compareId.
subdetId() ||
653 if (!itmCompare->recHit()->isValid())
continue;
674 previousTM = &(* itTrajMeas);
680 hits.push_back(hit->
clone());
688 else if(verdict==-2) hits.push_back(hit->
clone());
689 else if(verdict==-1){
691 hits.push_back(hit->
clone());
696 std::reverse(hits.begin(),hits.end());
709 for (std::vector<Rule>::const_iterator itr =
rules_.begin(), edr =
rules_.end(); itr != edr; ++itr) {
710 itr->apply(detid, verdict);
722 if(!
checkStoN(iSetup, detid, hit))hitresult=-5;
735 bool keepthishit=
true;
740 int subdet_cnt=
id.subdetId();
749 const std::type_info &
type =
typeid(*therechit);
753 if (hit!=0) cluster = &*(hit->
cluster());
755 edm::LogError(
"TrackerTrackHitFilter")<<
"TrackerTrackHitFilter::checkStoN : Unknown valid tracker hit in subdet " <<
id.subdetId()<<
"(detID="<<
id.rawId()<<
")\n ";
761 if (hit!=0) cluster = &*(hit->
cluster());
763 edm::LogError(
"TrackerTrackHitFilter")<<
"TrackerTrackHitFilter::checkStoN : Unknown valid tracker hit in subdet " <<
id.subdetId()<<
"(detID="<<
id.rawId()<<
")\n ";
771 throw cms::Exception(
"Unknown RecHit Type") <<
"RecHit of type " << type.name() <<
772 " not supported. (use c++filt to demangle the name)";
785 else if (subdet_cnt<=2){
804 int qbin =pixelhit->
qBin();
815 else {
std::cout<<
"HIT IN PIXEL ("<<subdet_cnt <<
") but PixelRecHit is EMPTY!!!"<<std::endl;}
831 bool corrcharge_ok=
true;
846 float angle=TMath::ASin(TMath::Abs(mom_z) /
sqrt(
pow(mom_x,2)+
pow(mom_y,2)+
pow(mom_z,2) ) );
854 if(hitpointer->isValid()){
858 float clust_alpha= atan2( mom_z, mom_x );
859 float clust_beta= atan2( mom_z, mom_y );
864 float clust_charge=pixelhit->
cluster()->charge();
865 float corr_clust_charge=clust_charge *
sqrt( 1.0 / ( 1.0/
pow(
tan(clust_alpha), 2 ) +
866 1.0/
pow(
tan(clust_beta ), 2 ) +
880 edm::LogWarning(
"TrackerTrackHitFilter") <<
"TSOS not valid ! Impossible to calculate track angle.";
883 return angle_ok&&corrcharge_ok;
892 bool corrcharge_ok=
false;
895 if(!hitpointer->isValid())
return corrcharge_ok;
898 return corrcharge_ok;
906 float clust_alpha= atan2( mom_z, mom_x );
907 float clust_beta= atan2( mom_z, mom_y );
912 float clust_charge=pixelhit->
cluster()->charge();
913 float corr_clust_charge=clust_charge *
sqrt( 1.0 / ( 1.0/
pow(
tan(clust_alpha), 2 ) +
914 1.0/
pow(
tan(clust_beta ), 2 ) +
920 return corrcharge_ok;
930 return tobId.
layer();
934 return tobId.
disk() + (3*(tobId.
side()-1));
938 return tibId.
layer();
942 return tobId.
layer();
946 return tobId.
wheel() + (9*(tobId.
side()-1));
950 return tobId.
wheel() + (3*(tobId.
side()-1));
ClusterRef cluster() const
T getParameter(std::string const &) const
EventNumber_t event() const
virtual void produce(edm::Event &iEvent, const edm::EventSetup &iSetup)
bool hasFilledProb() const
float clusterProbability(unsigned int flags=0) const
const Vector & momentum() const
track momentum vector
unsigned int layer() const
layer id
int layerFromId(const DetId &id) const
std::vector< double > subdetStoNlowcut_
#define DEFINE_FWK_MODULE(type)
LocalVector localDirection() const
TrackerTrackHitFilter(const edm::ParameterSet &iConfig)
std::vector< TrackCandidate > TrackCandidateCollection
bool checkStoN(const edm::EventSetup &iSetup, const DetId &id, const TrackingRecHit *therechit)
ConstRecHitPointer recHit() const
bool stripFrontInvalidHits_
const math::XYZPoint & outerPosition() const
position of the outermost hit
bool stripBackInvalidHits_
unsigned int layer() const
layer id
void produceFromTrajectory(const edm::EventSetup &iSetup, const Trajectory *itt, std::vector< TrackingRecHit * > &hits)
unsigned int side() const
positive or negative id
void apply(DetId detid, bool &verdict) const
const math::XYZPoint & innerPosition() const
position of the innermost hit
uint32_t rawId() const
get the raw id
U second(std::pair< T, U > const &p)
bool checkPXLCorrClustCharge(const TrajectoryMeasurement &meas)
std::vector< uint32_t > detsToIgnore_
DataContainer const & measurements() const
double PXLcorrClusChargeCut_
float signalOverNoise() const
std::vector< Rule > rules_
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Rule(const std::string &str)
int layer(DetId detid) const
void produceFromTrack(const edm::EventSetup &iSetup, const Track *itt, std::vector< TrackingRecHit * > &hits)
tuple TrackerTrackHitFilter
bool replaceWithInactiveHits_
unsigned int outerDetId() const
DetId of the detector on which surface the outermost state is located.
Tan< T >::type tan(const T &t)
ClusterRef cluster() const
TrajectoryStateOnSurface updatedState() const
trackingRecHit_iterator recHitsBegin() const
Iterator to first hit on the track.
int subdetId() const
get the contents of the subdetector field (not cast into any detector's numbering enum) ...
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
unsigned int disk() const
disk id
std::vector< double > subdetStoNhighcut_
virtual TrackingRecHit * clone() const =0
Detector identifier class for the strip tracker.
edm::RefToBase< TrajectorySeed > seedRef() const
edm::ESHandle< TransientTrackingRecHitBuilder > theBuilder
int checkHit(const edm::EventSetup &iSetup, const DetId &detid, const TrackingRecHit *hit)
edm::ESHandle< TrackerGeometry > theGeometry
unsigned int side() const
positive or negative id
ClusterRef cluster() const
key_type key() const
Accessor for product key.
void parseStoN(const std::string &str)
unsigned int wheel() const
wheel id
unsigned int layer() const
layer id
std::vector< int32_t > pxlTPLqBin_
bool checkHitAngle(const TrajectoryMeasurement &meas)
unsigned int side() const
positive or negative id
std::vector< bool > subdetStoN_
PropagationDirection seedDirection() const
direction of how the hits were sorted in the original seed
bool stripAllInvalidHits_
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.
DetId geographicalId() const
Detector det() const
get the detector field from this detid
TrackCandidate makeCandidate(const reco::Track &tk, std::vector< TrackingRecHit * >::iterator hitsBegin, std::vector< TrackingRecHit * >::iterator hitsEnd)
edm::ESHandle< MagneticField > theMagField
unsigned int innerDetId() const
DetId of the detector on which surface the innermost state is located.
std::string CMNSubtractionMode_
Power< A, B >::type pow(const A &a, const B &b)
unsigned int wheel() const
wheel id
T angle(T x1, T y1, T z1, T x2, T y2, T z2)
trackingRecHit_iterator recHitsEnd() const
Iterator to last hit on the track.