Skip to content

Commit

Permalink
- Configuration for run test
Browse files Browse the repository at this point in the history 
- scincillator configuration are added
  • Loading branch information
@quark2
quark2 committed Jul 31, 2017
1 parent aeb3119 commit 095cff4
Show file tree
Hide file tree
Showing 8 changed files with 325 additions and 41 deletions.
29 changes: 25 additions & 4 deletions EventFilter/GEMRawToDigi/plugins/GEMCosmicStandUnpacker.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -87,7 +87,9 @@ class GEMCosmicStandUnpacker : public edm::EDProducer {
std::vector<int> columnVector_;
std::vector<int> layerVector_;
std::vector<int> hitsVector_;
std::vector<int> nullstripVector_;
std::vector<unsigned long long> missingVector_;
std::vector<int> mulmissing_;
//std::vector<unsigned long long> missingPositionInVector_; // not useful


Expand Down Expand Up @@ -169,6 +171,7 @@ GEMCosmicStandUnpacker::beginRun(const edm::Run &run, const edm::EventSetup& iSe

for (unsigned int i=0; i<vfatVector_.size(); i++){
hitsVector_.push_back(0);
nullstripVector_.push_back(0);
}


Expand Down Expand Up @@ -345,13 +348,16 @@ GEMCosmicStandUnpacker::produce(edm::Event& iEvent, const edm::EventSetup& iSetu
if(crc!=crc_check) std::cout<<"DIFFERENT CRC :"<<crc<<" "<<crc_check<<std::endl;

bool Quality = (b1010==10) && (b1100==12) && (b1110==14) && (crc==crc_check) ;

if(!Quality && checkQualityEvent_) {std::cout << "Bad quality" << std::endl; continue;}

uint64_t converted=ChipID+0xf000;
bool foundChip=false;
int column=1;
int row=1;
int chamberPosition=1;
int slot=0;
int nNoVFAT = 0;
for(unsigned int i=0; i<vfatVector_.size(); i++) {
if( converted == vfatVector_[i]) {
foundChip=true;
Expand All @@ -360,15 +366,21 @@ GEMCosmicStandUnpacker::produce(edm::Event& iEvent, const edm::EventSetup& iSetu
row=rowVector_[i];
chamberPosition=layerVector_[i];
hitsVector_[i]++;
nNoVFAT = i;
}}
int schamberPosition=1+2*(row-1)+10*(column-1);
if (!foundChip) {
bool alreadyin=false;
for (unsigned int i=0; i<missingVector_.size(); i++) {
if(converted==missingVector_[i]) alreadyin=true;
if(converted==missingVector_[i]) {
alreadyin=true;
mulmissing_[ i ]++;
std::cout<<"##### Multiple VFAT not in the configuration!!!"<<std::endl;
}
}
if(!alreadyin) {
missingVector_.push_back(converted);
mulmissing_.push_back(1);
//missingPositionInVector_.push_back(k);
std::cout<<"Unpacked VFAT not in the configuration - double check the settings"<<std::endl;
std::cout<<" ---> VFAT--->"<<converted<<" (will only give this warning once)"<<std::endl;
Expand All @@ -382,10 +394,12 @@ GEMCosmicStandUnpacker::produce(edm::Event& iEvent, const edm::EventSetup& iSetu
std::cout<<std::dec<<" --> COLUMN = "<<column<<" ROW = "<<row<<" Layer:"<<chamberPosition<<" --> SC:"<<schamberPosition<<std::endl;
}

if(!Quality && checkQualityEvent_) continue;
//if(!Quality && checkQualityEvent_) continue;
if(!Quality && checkQualityEvent_) {std::cout << "Bad quality" << std::endl; continue;}

int bx=0;
uint8_t chan0xf = 0;
int nIsFilled = 0;

for(int chan = 0; chan < 128; ++chan) {

Expand All @@ -396,6 +410,7 @@ GEMCosmicStandUnpacker::produce(edm::Event& iEvent, const edm::EventSetup& iSetu
}

if(chan0xf==0) continue;
nIsFilled = 1;

GEMROmap::eCoord ec;
ec.chamberId=31;
Expand All @@ -414,6 +429,7 @@ GEMCosmicStandUnpacker::produce(edm::Event& iEvent, const edm::EventSetup& iSetu

int etaP=dc.etaId;

//std::cout << "etaId, chamberPos, schamberPos, etaId : " << etaP << ", " << chamberPosition << ", " << schamberPosition << ", " << etaP << std::endl;
if(etaP == 0) {
if(foundChip){
std::cout<<"WARNING: wrong digi! DoubleCheck the configuration"<<std::endl;
Expand All @@ -434,6 +450,11 @@ GEMCosmicStandUnpacker::produce(edm::Event& iEvent, const edm::EventSetup& iSetu
}
}
delete m_vfatdata;

if ( nIsFilled == 0 ) {
nullstripVector_[ nNoVFAT ]++;
//std::cout << "null strip" << std::endl;
}

}
std::fread(&m_word, sizeof(uint64_t), 1, m_file);
Expand Down Expand Up @@ -557,11 +578,11 @@ GEMCosmicStandUnpacker::endRun(edm::Run const&, edm::EventSetup const&)
std::cout<<"FILLED VFATS"<<std::endl;
std::cout<<"SLOT - VFAT - LAYER - COLUMN - ROW ---> Number of events with hits"<<std::endl;
for (unsigned int i=0; i<vfatVector_.size(); i++){
std::cout<<std::dec<<"\t"<<slotVector_.at(i)<<"\t"<<std::hex<<vfatVector_.at(i)<<"\t"<<layerVector_.at(i)<<"\t"<<columnVector_.at(i)<<"\t"<<rowVector_.at(i)<<" -----> "<<std::dec<<hitsVector_.at(i)<<std::endl;
std::cout<<std::dec<<"\t"<<slotVector_.at(i)<<"\t"<<std::hex<<vfatVector_.at(i)<<"\t"<<layerVector_.at(i)<<"\t"<<columnVector_.at(i)<<"\t"<<rowVector_.at(i)<<" -----> "<<std::dec<<hitsVector_.at(i)<<"\t("<<hitsVector_.at(i) - nullstripVector_.at(i)<<")"<<std::endl;
}
std::cout<<"MISSING VFATS"<<std::endl;
for (unsigned int i=0; i<missingVector_.size(); i++){
std::cout<<std::hex<<missingVector_[i]<<std::dec<<std::endl;
std::cout<<std::hex<<missingVector_[i]<<std::dec<<"\t("<<mulmissing_[ i ]<<")"<<std::endl;
}


Expand Down
8 changes: 7 additions & 1 deletion IOMC/ParticleGuns/src/FlatRandomPtGunProducer.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -62,7 +62,10 @@ void FlatRandomPtGunProducer::produce(Event &e, const EventSetup& es)
//
// 1st, primary vertex
//
HepMC::GenVertex* Vtx = new HepMC::GenVertex(HepMC::FourVector(0.,0.,0.));
//HepMC::GenVertex* Vtx = new HepMC::GenVertex(HepMC::FourVector(0.,0.,0.));
double dVx = CLHEP::RandFlat::shoot(engine, -1000.0, 1000.0) ;
double dVz = CLHEP::RandFlat::shoot(engine, -1000.0, 1000.0) ;
HepMC::GenVertex* Vtx = new HepMC::GenVertex(HepMC::FourVector(dVx,-100.0,dVz));

// loop over particles
//
Expand All @@ -81,6 +84,9 @@ void FlatRandomPtGunProducer::produce(Event &e, const EventSetup& es)
double px = pt*cos(phi) ;
double py = pt*sin(phi) ;
double pz = mom*cos(theta) ;
px = 0;
py = pt;
pz = 0;
double energy2= mom*mom + mass*mass ;
double energy = sqrt(energy2) ;
HepMC::FourVector p(px,py,pz,energy) ;
Expand Down
5 changes: 5 additions & 0 deletions Validation/GEMCR/interface/gemcrValidation.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -83,6 +83,11 @@ class gemcrValidation : public GEMBaseValidation
std::auto_ptr<std::vector<TrajectorySeed> > findSeeds(MuonTransientTrackingRecHit::MuonRecHitContainer &muRecHits);
Trajectory makeTrajectory(TrajectorySeed seed, MuonTransientTrackingRecHit::MuonRecHitContainer &muRecHits, std::vector<GEMChamber> gemChambers, GEMChamber testChamber);
edm::EDGetToken InputTagToken_, InputTagToken_RH, InputTagToken_TR, InputTagToken_TS, InputTagToken_DG;

float fScinHPosY, fScinHLeft, fScinHRight, fScinHTop, fScinHBottom;
float fScinLPosY, fScinLLeft, fScinLRight, fScinLTop, fScinLBottom;

bool isPassedScincillators(GlobalPoint trajGP1, GlobalPoint trajGP2);
};

#endif
64 changes: 53 additions & 11 deletions Validation/GEMCR/src/gemcrValidation.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -31,13 +31,6 @@

#include "DataFormats/GEMDigi/interface/GEMDigiCollection.h"

#include "FWCore/ServiceRegistry/interface/Service.h"
#include "FWCore/Utilities/interface/RandomNumberGenerator.h"
#include "CLHEP/Random/RandFlat.h"




#include <iomanip>

#include <TCanvas.h>
Expand All @@ -58,6 +51,19 @@ gemcrValidation::gemcrValidation(const edm::ParameterSet& cfg): GEMBaseValidatio
trackChi2 = cfg.getParameter<double>("trackChi2");
trackResY = cfg.getParameter<double>("trackResY");
trackResX = cfg.getParameter<double>("trackResX");

fScinHPosY = cfg.getParameter<double>("ScincilUpperY");
fScinHLeft = cfg.getParameter<double>("ScincilUpperLeft");
fScinHRight = cfg.getParameter<double>("ScincilUpperRight");
fScinHTop = cfg.getParameter<double>("ScincilUpperTop");
fScinHBottom = cfg.getParameter<double>("ScincilUpperBottom");

fScinLPosY = cfg.getParameter<double>("ScincilLowerY");
fScinLLeft = cfg.getParameter<double>("ScincilLowerLeft");
fScinLRight = cfg.getParameter<double>("ScincilLowerRight");
fScinLTop = cfg.getParameter<double>("ScincilLowerTop");
fScinLBottom = cfg.getParameter<double>("ScincilLowerBottom");

edm::ParameterSet smootherPSet = cfg.getParameter<edm::ParameterSet>("MuonSmootherParameters");
theSmoother = new CosmicMuonSmoother(smootherPSet, theService);
theUpdator = new KFUpdator();
Expand Down Expand Up @@ -91,7 +97,7 @@ void gemcrValidation::bookHistograms(DQMStore::IBooker & ibooker, edm::Run const
firedMul = ibooker.book1D("firedMul","fired chamber multiplicity",n_ch+1,0,n_ch+1);
firedChamber = ibooker.book1D("firedChamber", "fired chamber",n_ch,0,n_ch);

tr_chamber = ibooker.book1D("tr_eff_ch", "tr rec /chamber",n_ch,0,n_ch);
tr_chamber = ibooker.book1D("tr_eff_ch", "tr rec /chamber",n_ch,1,n_ch);
th_chamber = ibooker.book1D("th_eff_ch", "tr hit/chamber",n_ch,0,n_ch);
rh_chamber = ibooker.book1D("rh_eff_ch", "rec hit/chamber",n_ch,0,n_ch);
rh1_chamber = ibooker.book1D("rh1_chamber", "all recHits",n_ch,0,n_ch);
Expand Down Expand Up @@ -251,26 +257,55 @@ Trajectory gemcrValidation::makeTrajectory(TrajectorySeed seed, MuonTransientTra
return fitted.front();
}

void gemcrValidation::analyze(const edm::Event& e, const edm::EventSetup& iSetup){

bool gemcrValidation::isPassedScincillators(GlobalPoint trajGP1, GlobalPoint trajGP2) {
// Getting the direction of trajectory by using the hits in the seed
float fTrajDirX = trajGP2.x() - trajGP1.x();
float fTrajDirY = trajGP2.y() - trajGP1.y();
float fTrajDirZ = trajGP2.z() - trajGP1.z();

// Time to solve 1-order equation
float fTHTraj = ( fScinHPosY - trajGP1.y() ) / fTrajDirY;
float fTLTraj = ( fScinLPosY - trajGP1.y() ) / fTrajDirY;

// Finding the coordinates of the hit position on the y-plane containing the upper scincillator
float fXHitH = trajGP1.x() + fTHTraj * fTrajDirX;
float fZHitH = trajGP1.z() + fTHTraj * fTrajDirZ;

// Checking whether the hit position in the upper scincillator
if ( !( ( fScinHLeft <= fXHitH && fXHitH <= fScinHRight ) && ( fScinHTop <= fZHitH && fZHitH <= fScinHBottom ) ) ) return false;

// Finding the coordinates of the hit position on the y-plane containing the lower scincillator
float fXHitL = trajGP1.x() + fTLTraj * fTrajDirX;
float fZHitL = trajGP1.z() + fTLTraj * fTrajDirZ;

// Checking whether the hit position in the lower scincillator
if ( !( ( fScinLLeft <= fXHitL && fXHitL <= fScinLRight ) && ( fScinLTop <= fZHitL && fZHitL <= fScinLBottom ) ) ) return false;

return true;
}

edm::Service<edm::RandomNumberGenerator> rng;
CLHEP::HepRandomEngine* engine = &rng->getEngine(e.streamID());

void gemcrValidation::analyze(const edm::Event& e, const edm::EventSetup& iSetup){

//std::cout << "analyze on!" << std::endl;
theService->update(iSetup);

edm::Handle<GEMRecHitCollection> gemRecHits;
if (!isMC){
edm::Handle<GEMDigiCollection> digis;
e.getByToken( this->InputTagToken_DG, digis);
int nNumDigi = 0;
for (GEMDigiCollection::DigiRangeIterator gemdgIt = digis->begin(); gemdgIt != digis->end(); ++gemdgIt){
nNumDigi++;
const GEMDetId& gemId = (*gemdgIt).first;
int index = findIndex(gemId);
const GEMDigiCollection::Range& range = (*gemdgIt).second;
for ( auto digi = range.first; digi != range.second; ++digi ) {
gem_chamber_digi_digi[index]->Fill(digi->strip(),gemId.roll());
}
}
//std::cout << "num of digi : " << nNumDigi << std::endl;
}
e.getByToken( this->InputTagToken_RH, gemRecHits);
if (!gemRecHits.isValid()) {
Expand All @@ -285,6 +320,7 @@ void gemcrValidation::analyze(const edm::Event& e, const edm::EventSetup& iSetup
firedCh.push_back(0);
rMul.push_back(0);
}
//std::cout << "num of hit : " << gemRecHits->size() << std::endl;
for (GEMRecHitCollection::const_iterator recHit = gemRecHits->begin(); recHit != gemRecHits->end(); ++recHit){

Float_t rh_l_x = recHit->localPosition().x();
Expand All @@ -299,6 +335,7 @@ void gemcrValidation::analyze(const edm::Event& e, const edm::EventSetup& iSetup
//checkRH[index] = 1;
Short_t rh_roll = (Short_t) id.roll();
LocalPoint recHitLP = recHit->localPosition();
//printf("%i, %lf, %i, %i\n", index, rh_l_x, clusterSize, firstClusterStrip);

if ( GEMGeometry_->idToDet((*recHit).gemId()) == nullptr) {
std::cout<<"This gem recHit did not matched with GEMGeometry."<<std::endl;
Expand Down Expand Up @@ -433,10 +470,15 @@ void gemcrValidation::analyze(const edm::Event& e, const edm::EventSetup& iSetup
if (!bestTrajectory.isValid()) continue; //{cout<<"no Best Trajectory" << endl; continue;}
gem_chamber_track[findIndex(tch.id())]->Fill(1.5);
if (maxChi2 > trackChi2) continue;

// Checking whether the trajectory passed scincillators
if ( !isPassedScincillators(bestSeed.recHits().first->globalPosition(), bestSeed.recHits().second->globalPosition()) ) continue;

trajectoryh->Fill(3,1);
gem_chamber_bestChi2[findIndex(tch.id())]->Fill(maxChi2);
//cout << maxChi2 << endl;
gem_chamber_track[findIndex(tch.id())]->Fill(2.5);

PTrajectoryStateOnDet ptsd1(bestSeed.startingState());
DetId did(ptsd1.detId());
const BoundPlane& bp = theService->trackingGeometry()->idToDet(did)->surface();
Expand Down
Loading

0 comments on commit 095cff4

Please sign in to comment.