Skip to content

Commit

Permalink
Merge pull request #92 from cms-sw/CMSSW_6_2_X_SLHC
Browse files Browse the repository at this point in the history
Cmssw 6 2 x slhc
  • Loading branch information
jshlee committed Mar 4, 2014
2 parents fe65259 + c81720e commit 50133bb
Show file tree
Hide file tree
Showing 4 changed files with 115 additions and 158 deletions.
102 changes: 55 additions & 47 deletions L1Trigger/TrackTrigger/plugins/TTStubBuilder.h
Original file line number Diff line number Diff line change
@@ -1,17 +1,17 @@
/*! \class TTStubBuilder
* \brief Plugin to load the Stub finding algorithm and produce the
* collection of Stubs that goes in the event content.
* \details After moving from SimDataFormats to DataFormats,
* the template structure of the class was maintained
* in order to accomodate any types other than PixelDigis
* in case there is such a need in the future.
*
* \author Andrew W. Rose
* \author Nicola Pozzobon
* \author Ivan Reid
* \date 2013, Jul 18
*
*/
/*! \class TTStubBuilder
* \brief Plugin to load the Stub finding algorithm and produce the
* collection of Stubs that goes in the event content.
* \details After moving from SimDataFormats to DataFormats,
* the template structure of the class was maintained
* in order to accomodate any types other than PixelDigis
* in case there is such a need in the future.
*
* \author Andrew W. Rose
* \author Nicola Pozzobon
* \author Ivan Reid
* \date 2013, Jul 18
*
*/

#ifndef L1_TRACK_TRIGGER_STUB_BUILDER_H
#define L1_TRACK_TRIGGER_STUB_BUILDER_H
Expand Down Expand Up @@ -48,9 +48,9 @@ class TTStubBuilder : public edm::EDProducer

private:
/// Data members
const StackedTrackerGeometry *theStackedTracker;
const StackedTrackerGeometry *theStackedTracker;
edm::ESHandle< TTStubAlgorithm< T > > theStubFindingAlgoHandle;
edm::InputTag TTClustersInputTag;
edm::InputTag TTClustersInputTag;

/// Mandatory methods
virtual void beginRun( const edm::Run& run, const edm::EventSetup& iSetup );
Expand All @@ -63,12 +63,12 @@ class TTStubBuilder : public edm::EDProducer

}; /// Close class

/*! \brief Implementation of methods
* \details Here, in the header file, the methods which do not depend
* on the specific type <T> that can fit the template.
* Other methods, with type-specific features, are implemented
* in the source file.
*/
/*! \brief Implementation of methods
* \details Here, in the header file, the methods which do not depend
* on the specific type <T> that can fit the template.
* Other methods, with type-specific features, are implemented
* in the source file.
*/

/// Constructors
template< typename T >
Expand All @@ -93,7 +93,7 @@ void TTStubBuilder< T >::beginRun( const edm::Run& run, const edm::EventSetup& i
iSetup.get< StackedTrackerGeometryRecord >().get( StackedTrackerGeomHandle );
theStackedTracker = StackedTrackerGeomHandle.product();

/// Get the stub finding algorithm
/// Get the stub finding algorithm
iSetup.get< TTStubAlgorithmRecord >().get( theStubFindingAlgoHandle );

/// Print some information when loaded
Expand All @@ -111,7 +111,7 @@ void TTStubBuilder< T >::endRun( const edm::Run& run, const edm::EventSetup& iSe
/// Implement the producer
template< typename T >
void TTStubBuilder< T >::produce( edm::Event& iEvent, const edm::EventSetup& iSetup )
{
{
/// Prepare output
std::auto_ptr< edmNew::DetSetVector< TTCluster< T > > > TTClusterDSVForOutput( new edmNew::DetSetVector< TTCluster< T > > );
std::auto_ptr< edmNew::DetSetVector< TTStub< T > > > TTStubDSVForOutputTemp( new edmNew::DetSetVector< TTStub< T > > );
Expand Down Expand Up @@ -186,7 +186,7 @@ void TTStubBuilder< T >::produce( edm::Event& iEvent, const edm::EventSetup& iSe
const GeomDetUnit* det0 = theStackedTracker->idToDetUnit( Id, 0 );
const PixelGeomDetUnit* pix0 = dynamic_cast< const PixelGeomDetUnit* >( det0 );
const PixelTopology* top0 = dynamic_cast< const PixelTopology* >( &(pix0->specificTopology()) );
const int chipSize = 2 * top0->rowsperroc(); /// Need to find ASIC size in half-strip units
const int chipSize = 2 * top0->rowsperroc(); /// Need to find ASIC size in half-strip units
std::map< int, std::vector< TTStub< T > > > moduleStubs; /// Temporary storage for stubs before max check

/// Loop over pairs of Clusters
Expand All @@ -206,7 +206,7 @@ void TTStubBuilder< T >::produce( edm::Event& iEvent, const edm::EventSetup& iSe
/// Check for compatibility
bool thisConfirmation = false;
int thisDisplacement = 999999;
int thisOffset = 0;
int thisOffset = 0;

theStubFindingAlgoHandle->PatternHitCorrelation( thisConfirmation, thisDisplacement, thisOffset, tempTTStub );

Expand All @@ -229,7 +229,7 @@ void TTStubBuilder< T >::produce( edm::Event& iEvent, const edm::EventSetup& iSe
/// This means that only some of them do
/// Put in the temporary output
int chip = tempTTStub.getTriggerPosition() / chipSize; /// Find out which ASIC
if ( moduleStubs.find( chip ) == moduleStubs.end() ) /// Already a stub for this ASIC?
if ( moduleStubs.find( chip ) == moduleStubs.end() ) /// Already a stub for this ASIC?
{
/// No, so new entry
std::vector< TTStub< T > > tempStubs;
Expand All @@ -244,10 +244,13 @@ void TTStubBuilder< T >::produce( edm::Event& iEvent, const edm::EventSetup& iSe
}
}
} /// Stub accepted
else
{
tempRejected->push_back( tempTTStub );
} /// Stub rejected
/* NP 2014 02 25
* this is commented to avoid memory exhaustion in hi PU events
else
{
tempRejected->push_back( tempTTStub );
} /// Stub rejected
*/
} /// End of nested loop
} /// End of loop over pairs of Clusters

Expand Down Expand Up @@ -286,11 +289,14 @@ void TTStubBuilder< T >::produce( edm::Event& iEvent, const edm::EventSetup& iSe
tempOuter->push_back( *(is.second[bendMap[i].first].getClusterRef(1)) );
tempOutput->push_back( is.second[bendMap[i].first] );
}
for ( unsigned int i = maxStubs; i < is.second.size(); ++i )
{
/// Reject the rest
tempRejected->push_back( is.second[bendMap[i].first] );
}
/* NP 2014 02 25
* this is commented to avoid memory exhaustion in hi PU events
for ( unsigned int i = maxStubs; i < is.second.size(); ++i )
{
/// Reject the rest
tempRejected->push_back( is.second[bendMap[i].first] );
}
*/
}
} /// End of loop over temp output
} /// End store only the selected stubs if max no. stub/ROC is set
Expand Down Expand Up @@ -329,16 +335,19 @@ void TTStubBuilder< T >::produce( edm::Event& iEvent, const edm::EventSetup& iSe
tempOutputFiller.abort();
}

if ( tempRejected->size() > 0 )
{
typename edmNew::DetSetVector< TTStub< T > >::FastFiller rejectedOutputFiller( *TTStubDSVForOutputRejected, DetId(Id.rawId()) );
for ( unsigned int m = 0; m < tempRejected->size(); m++ )
{
rejectedOutputFiller.push_back( tempRejected->at(m) );
}
if ( rejectedOutputFiller.empty() )
rejectedOutputFiller.abort();
}
/* NP 2014 02 25
* this is commented to avoid memory exhaustion in hi PU events
if ( tempRejected->size() > 0 )
{
typename edmNew::DetSetVector< TTStub< T > >::FastFiller rejectedOutputFiller( *TTStubDSVForOutputRejected, DetId(Id.rawId()) );
for ( unsigned int m = 0; m < tempRejected->size(); m++ )
{
rejectedOutputFiller.push_back( tempRejected->at(m) );
}
if ( rejectedOutputFiller.empty() )
rejectedOutputFiller.abort();
}
*/

} /// End of loop over detector elements

Expand Down Expand Up @@ -461,4 +470,3 @@ bool TTStubBuilder< T >::SortStubBendPairs( const std::pair< unsigned int, doubl
}

#endif

Original file line number Diff line number Diff line change
Expand Up @@ -57,10 +57,6 @@ using namespace reco;
using namespace l1extra;


//bool myfunction (float i,float j) { return (i>j); }

//bool sortTLorentz (TLorentzVector i,TLorentzVector j) { return ( i.Pt()>j.Pt() ); }


class L1CalibFilterTowerJetProducer : public edm::EDProducer {
public:
Expand All @@ -79,24 +75,11 @@ class L1CalibFilterTowerJetProducer : public edm::EDProducer {
virtual void beginLuminosityBlock(edm::LuminosityBlock&, edm::EventSetup const&);
virtual void endLuminosityBlock(edm::LuminosityBlock&, edm::EventSetup const&);

float get_rho(double L1rho);

//fwd calibration: V ROUGH (only to L1extra particles)
// double rough_ptcal(double pt);

double Median(vector<double> aVec);

void TMVA_calibration();
// ----------member data ---------------------------
ParameterSet conf_;

ifstream indata;


edm::FileInPath inMVAweights_edm;
float l1Pt, l1Eta;
TMVA::Reader *reader;
float val_pt_cal(float l1pt, float l1eta);
};


Expand All @@ -111,8 +94,6 @@ conf_(iConfig)
produces< L1EtMissParticleCollection >( "TowerMHT" ) ;
produces<double>("TowerHT");

//look up tables
inMVAweights_edm = iConfig.getParameter<edm::FileInPath> ("inMVA_weights_file");

}

Expand Down Expand Up @@ -174,8 +155,9 @@ L1CalibFilterTowerJetProducer::produce(edm::Event& iEvent, const edm::EventSetup
float l1wEta_ = il1->WeightedEta();
float l1wPhi_ = il1->WeightedPhi() ;

//Get the calibration factors from the MVA lookup table
float cal_Pt_ = val_pt_cal( l1Pt_ , l1wEta_) * l1Pt_;

// Currently no calibration is applied
double cal_Pt_ = l1Pt_;

math::PtEtaPhiMLorentzVector p4;

Expand All @@ -190,22 +172,16 @@ L1CalibFilterTowerJetProducer::produce(edm::Event& iEvent, const edm::EventSetup
if( cal_Pt_>15 ) mht+=upgrade_jet;

// add jet to L1Extra list
outputExtraCen->push_back( L1JetParticle( math::PtEtaPhiMLorentzVector( cal_Pt_,
l1wEta_,
l1wPhi_,
0. ),
outputExtraCen->push_back( L1JetParticle( math::PtEtaPhiMLorentzVector( cal_Pt_, l1wEta_, l1wPhi_, 0. ),
Ref< L1GctJetCandCollection >(), 0 )
);
}




// create L1Extra object
math::PtEtaPhiMLorentzVector p4tmp = math::PtEtaPhiMLorentzVector( mht.pt(),
0.,
mht.phi(),
0. ) ;
// create L1Extra object
math::PtEtaPhiMLorentzVector p4tmp = math::PtEtaPhiMLorentzVector( mht.pt(), 0., mht.phi(), 0. ) ;

L1EtMissParticle l1extraMHT(p4tmp,
L1EtMissParticle::kMHT,
Expand Down Expand Up @@ -233,9 +209,7 @@ L1CalibFilterTowerJetProducer::produce(edm::Event& iEvent, const edm::EventSetup

// ------------ method called once each job just before starting event loop ------------
void
L1CalibFilterTowerJetProducer::beginJob()
{

L1CalibFilterTowerJetProducer::beginJob() {
}

// ------------ method called once each job just after ending the event loop ------------
Expand All @@ -247,9 +221,6 @@ L1CalibFilterTowerJetProducer::endJob() {
void
L1CalibFilterTowerJetProducer::beginRun(edm::Run&, edm::EventSetup const&)
{

//read in calibration for pt from TMVA
TMVA_calibration();
}

// ------------ method called when ending the processing of a run ------------
Expand Down Expand Up @@ -280,40 +251,8 @@ L1CalibFilterTowerJetProducer::fillDescriptions(edm::ConfigurationDescriptions&
descriptions.addDefault(desc);
}

//
// member functions
//

void L1CalibFilterTowerJetProducer::TMVA_calibration()
{
cout<<"Getting lookup from MVA"<<endl;

reader = new TMVA::Reader("!Color:Silent");
reader->AddVariable( "l1Pt", &l1Pt);
reader->AddVariable( "l1Eta", &l1Eta);
// reader->AddVariable( "l1Phi", &l1Phi);
// reader->AddVariable( "MVA_Rho", &MVA_Rho );

cout<<"Booking MVA reader"<<endl;

reader->BookMVA("BDT method",inMVAweights_edm.fullPath().c_str());

cout<<"MVA reader booked: start processing events."<<endl;

}

float
L1CalibFilterTowerJetProducer::val_pt_cal(float l1pt, float l1eta)
{
l1Eta=l1eta;
l1Pt=l1pt;

Float_t val = (reader->EvaluateRegression(TString("BDT method") ))[0];

//cout<<"l1 pt: " << l1pt <<" corr_pt_1 "<<corr_pt_1<<endl;
return val;

}



Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -112,9 +112,6 @@ conf_(iConfig)
{

produces<double>("Rho");

//look up tables
inRhoData_edm = iConfig.getParameter<edm::FileInPath> ("inRhodata_file");

}

Expand Down Expand Up @@ -169,12 +166,10 @@ L1TowerJetPUEstimator::produce(edm::Event& iEvent, const edm::EventSetup& iSetup
//apply calibration to the raw rho: should we do this at this stage?
//not sure of the effect on high PU data

double cal_rhoL1 = raw_rho2 * get_rho(raw_rho2);

///////////////////////////////////////////////////
// SET VALUE OF RHO
///////////////////////////////////////////////////
outrho=cal_rhoL1;
outrho=raw_rho2;

*outRho = outrho;

Expand All @@ -200,31 +195,6 @@ L1TowerJetPUEstimator::endJob() {
void
L1TowerJetPUEstimator::beginRun(edm::Run&, edm::EventSetup const&)
{

//read in calibration for rho lookup table

inrhodata.open(inRhoData_edm.fullPath().c_str());
if(!inrhodata) cerr<<" unable to open rho lookup file. "<<endl;


//read into a vector
pair<double, double> rho_cal;
double L1rho_(9999), calFac_(9999);
while ( !inrhodata.eof() ) { // keep reading until end-of-file
// sets EOF flag if no value found
inrhodata >> L1rho_ >> calFac_ ;

rho_cal.first = L1rho_;
rho_cal.second= calFac_;

rho_cal_vec.push_back(rho_cal);
}
inrhodata.close();

cout<<" Read in rho lookup table"<<endl;



}

// ------------ method called when ending the processing of a run ------------
Expand Down Expand Up @@ -277,17 +247,5 @@ double L1TowerJetPUEstimator::Median( vector<double> aVec){
}



double L1TowerJetPUEstimator::get_rho(double L1_rho)
{

//get the rho multiplication factor:
//L1_rho * 2 gets the array index
if(L1_rho<=40.5) return rho_cal_vec[L1_rho*2].second;
//for L1 rho > 40.5, calibration flattens out
else return 1.44576;

}

//define this as a plug-in
DEFINE_FWK_MODULE(L1TowerJetPUEstimator);
Loading

0 comments on commit 50133bb

Please sign in to comment.