Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Updated code merging request #246

Merged
merged 39 commits into from
Feb 15, 2018
Merged
Show file tree
Hide file tree
Changes from 26 commits
Commits
Show all changes
39 commits
Select commit Hold shift + click to select a range
2bd0de9
Merged GEMCosmicMuon from repository gem-sw
hyunyong Dec 13, 2016
b9ed997
Merge branch 'GEMCosmicMuon' of https://github.com/hyunyong/cmssw int…
hyunyong Dec 14, 2016
d5a1454
Merge branch 'GEMCosmicMuon' of https://github.com/hyunyong/cmssw int…
hyunyong Dec 14, 2016
bed60be
Merge branch 'GEMCosmicMuon' of https://github.com/gem-sw/cmssw into …
hyunyong Dec 17, 2016
29d3382
Merge branch 'GEMCosmicMuon' of https://github.com/hyunyong/cmssw int…
hyunyong Dec 18, 2016
ecbc708
Merge branch 'GEMCosmicMuon' of https://github.com/hyunyong/cmssw int…
hyunyong Dec 18, 2016
af425ec
update eff cal
hyunyong Jun 30, 2017
aeb3119
Merged GEMCosmicMuon from repository hyunyong
quark2 Jul 13, 2017
095cff4
- Configuration for run test
quark2 Jul 31, 2017
5daff45
some clean up
quark2 Aug 2, 2017
f310ffe
- Scintillator geometry is added on
quark2 Sep 9, 2017
9901e32
- Fixed something (names, etc.) of scintillators
quark2 Dec 3, 2017
e8eeaa7
- cos^2 distribution is applied
quark2 Dec 15, 2017
f176d25
The first version on top-bottom seeds
quark2 Dec 18, 2017
b3234c7
Completed almost developing the top-bottom seed method (one of import…
quark2 Jan 2, 2018
4d3750f
Completion of top-bottom seed method ; theUpdator->update() in makeTr…
quark2 Jan 2, 2018
9e26b6c
The part of reconstruction of trajectories is separated into RecoMuon…
quark2 Jan 3, 2018
ca150a1
cos^2 distrib, scintillator position modified, minor changes in *_sim.py
giovanni-mocellin Jan 3, 2018
ec3ad4d
Minor bug fixing
giovanni-mocellin Jan 3, 2018
52c0d9e
Removed some output files
giovanni-mocellin Jan 3, 2018
fcf24c2
Merge pull request #1 from giovanni-mocellin/CosmicStandModified
quark2 Jan 4, 2018
edb00a8
Getting track information
quark2 Jan 4, 2018
83085a0
Merge branch 'P5dataDQM' of https://github.com/quark2/cmssw into P5da…
quark2 Jan 4, 2018
cd54f95
Reduced drop of efficiencies from geometrical effect, some modificati…
quark2 Jan 10, 2018
685b0fe
Merged P5dataDQM from repository quark2 with cms-merge-topic
giovanni-mocellin Feb 6, 2018
40c63b4
Geometry modified, chi2 closest to 1
giovanni-mocellin Feb 8, 2018
bb644e2
Merged CosmicRayStandQC8 from repository giovanni-mocellin with cms-m…
kcc5837 Feb 8, 2018
10d244e
test
kcc5837 Feb 8, 2018
e1fb394
Merge pull request #4 from kcc5837/QC8analysis
giovanni-mocellin Feb 8, 2018
f81c47c
Centralized the gemVecChamType, added the script for the geometry cha…
giovanni-mocellin Feb 9, 2018
27228ae
Errors ironed out - still to be tested again
giovanni-mocellin Feb 12, 2018
f3749c3
Mods in sim.py, RecoMuon and Validation - exclude the chambers is not…
giovanni-mocellin Feb 12, 2018
8a71679
Added geometry for "removed" chambers: they are positioned outside th…
giovanni-mocellin Feb 12, 2018
e461393
Modified the position of removed chambers - AND this method is working!
giovanni-mocellin Feb 13, 2018
99ac22f
Plot gemcr fixed ( x -> -x ) otherwise the disposition of the chamber…
giovanni-mocellin Feb 13, 2018
80fbb93
Merged with Byeonghak´s code
giovanni-mocellin Feb 13, 2018
4fb0852
Forgot to add the definition of a histogram - Now added
giovanni-mocellin Feb 13, 2018
37ae924
Geometry modification: chambers moved in along-chamber-direction -> T…
giovanni-mocellin Feb 15, 2018
339b908
Merge branch 'GEMCosmicMuon' into CosmicRayStandQC8
giovanni-mocellin Feb 15, 2018
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
29 changes: 25 additions & 4 deletions EventFilter/GEMRawToDigi/plugins/GEMCosmicStandUnpacker.cc
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: 6 additions & 2 deletions IOMC/ParticleGuns/interface/BaseFlatGunProducer.h
Original file line number Diff line number Diff line change
Expand Up @@ -46,10 +46,14 @@ namespace edm

// gun particle(s) characteristics
std::vector<int> fPartIDs ;
double fMinEta ;
double fMaxEta ;
double fMinPhi ;
double fMaxPhi ;
double fMinEta ;
double fMaxEta ;
double fMinTheta ;
double fMaxTheta ;
bool fIsThetaFlat ; // If 'True': theta distribution is flat. If 'False': theta distribution is a cos^2


// the event format itself
HepMC::GenEvent* fEvt;
Expand Down
13 changes: 8 additions & 5 deletions IOMC/ParticleGuns/src/BaseFlatGunProducer.cc
Original file line number Diff line number Diff line change
Expand Up @@ -45,11 +45,14 @@ BaseFlatGunProducer::BaseFlatGunProducer( const ParameterSet& pset ) :
// it looks like it's NOT even necessary to check if it is,
// before trying to extract parameters - if it is empty,
// the default values seem to be taken
fPartIDs = pgun_params.getParameter< vector<int> >("PartID");
fMinEta = pgun_params.getParameter<double>("MinEta");
fMaxEta = pgun_params.getParameter<double>("MaxEta");
fMinPhi = pgun_params.getParameter<double>("MinPhi");
fMaxPhi = pgun_params.getParameter<double>("MaxPhi");
fPartIDs = pgun_params.getParameter< vector<int> >("PartID");
fMinPhi = pgun_params.getParameter<double>("MinPhi");
fMaxPhi = pgun_params.getParameter<double>("MaxPhi");
//fMinEta = pgun_params.getParameter<double>("MinEta");
//fMaxEta = pgun_params.getParameter<double>("MaxEta");
fMinTheta = pgun_params.getParameter<double>("MinTheta");
fMaxTheta = pgun_params.getParameter<double>("MaxTheta");
fIsThetaFlat = pgun_params.getParameter<bool>("IsThetaFlat"); // If 'True': theta distribution is flat. If 'False': theta distribution is a cos^2

//
//fPDGTablePath = "/afs/cern.ch/sw/lcg/external/clhep/1.9.2.1/slc3_ia32_gcc323/data/HepPDT/" ;
Expand Down
96 changes: 82 additions & 14 deletions IOMC/ParticleGuns/src/FlatRandomPtGunProducer.cc
Original file line number Diff line number Diff line change
Expand Up @@ -40,6 +40,34 @@ FlatRandomPtGunProducer::~FlatRandomPtGunProducer()
// no need to cleanup GenEvent memory - done in HepMCProduct
}


bool myIsMuonPassScint(double dVx, double dVy, double dVz, double dPx, double dPy, double dPz) {
// To test the drop-down of efficiency at edges, we can set the cut looser
double ScintilXMin = -1000.0;
double ScintilXMax = 1000.0;
double ScintilZMin = -605.6;
double ScintilZMax = 950.0;

double ScintilLowerY = -114.85;
double ScintilUpperY = 1540.15;

double dTLower = ( ScintilLowerY - dVy ) / dPy;
double dXLower = dVx + dTLower * dPx;
double dZLower = dVz + dTLower * dPz;

double dTUpper = ( ScintilUpperY - dVy ) / dPy;
double dXUpper = dVx + dTUpper * dPx;
double dZUpper = dVz + dTUpper * dPz;

if (( ScintilXMin <= dXLower && dXLower <= ScintilXMax && ScintilZMin <= dZLower && dZLower <= ScintilZMax ) &&
( ScintilXMin <= dXUpper && dXUpper <= ScintilXMax && ScintilZMin <= dZUpper && dZUpper <= ScintilZMax ))
{
return true;
}

else return false;
}

void FlatRandomPtGunProducer::produce(Event &e, const EventSetup& es)
{
edm::Service<edm::RandomNumberGenerator> rng;
Expand All @@ -62,30 +90,71 @@ 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, -1200.0, 1200.0) ;
double dVy = -100.0;
double dVz = CLHEP::RandFlat::shoot(engine, -650.0, 750.0) ;
HepMC::GenVertex* Vtx = new HepMC::GenVertex(HepMC::FourVector(dVx,dVy,dVz));*/
double dVx;
double dVy = 1540.15; // same Y as the upper scintillator
double dVz;
HepMC::GenVertex* Vtx = NULL;

// loop over particles
//
int barcode = 1 ;
for (unsigned int ip=0; ip<fPartIDs.size(); ++ip)
{
double pt = CLHEP::RandFlat::shoot(engine, fMinPt, fMaxPt) ;
double eta = CLHEP::RandFlat::shoot(engine, fMinEta, fMaxEta) ;
double phi = CLHEP::RandFlat::shoot(engine, fMinPhi, fMaxPhi) ;
double px, py, pz, mom;
double phi, theta;
int j = 0;

while (j < 10000) // j < 10000 to avoid too long computational time
{

dVx = CLHEP::RandFlat::shoot(engine, -1000.0, 1000.0) ;
dVz = CLHEP::RandFlat::shoot(engine, -605.6, 950.0) ;

mom = CLHEP::RandFlat::shoot(engine, fMinPt, fMaxPt) ;
phi = CLHEP::RandFlat::shoot(engine, fMinPhi, fMaxPhi) ;
theta = 0;

if (fIsThetaFlat)
{
theta = CLHEP::RandFlat::shoot(engine, fMinTheta, fMaxTheta);
}

if (!fIsThetaFlat)
{
double u = CLHEP::RandFlat::shoot(engine, 0.0, 0.785398); // u = Uniform[0;Pi/4]
theta = 0;
while(abs(u-(0.5*theta+0.25*sin(2*theta)))>0.000015)
{
theta+=0.00001;
}
}

px = mom*sin(theta)*cos(phi) ;
pz = mom*sin(theta)*sin(phi) ;
py = -mom*cos(theta) ; // with the - sign, the muons are going downwards: falling from the sky


if ( myIsMuonPassScint(dVx, dVy, dVz, px, py, pz) == true ) break; // muon passing through both the scintillators => valid: the loop can be stopped

j++;

}

int PartID = fPartIDs[ip] ;
const HepPDT::ParticleData*
PData = fPDGTable->particle(HepPDT::ParticleID(abs(PartID))) ;
double mass = PData->mass().value() ;
double theta = 2.*atan(exp(-eta)) ;
double mom = pt/sin(theta) ;
double px = pt*cos(phi) ;
double py = pt*sin(phi) ;
double pz = mom*cos(theta) ;
Vtx = new HepMC::GenVertex(HepMC::FourVector(dVx,dVy,dVz));

double energy2= mom*mom + mass*mass ;
double energy = sqrt(energy2) ;
double energy = sqrt(energy2) ;
HepMC::FourVector p(px,py,pz,energy) ;
HepMC::GenParticle* Part =
new HepMC::GenParticle(p,PartID,1);
HepMC::GenParticle* Part = new HepMC::GenParticle(p,PartID,1);
Part->suggest_barcode( barcode ) ;
barcode++ ;
Vtx->add_particle_out(Part);
Expand All @@ -98,8 +167,7 @@ void FlatRandomPtGunProducer::produce(Event &e, const EventSetup& es)
{
APartID = PartID ;
}
HepMC::GenParticle* APart =
new HepMC::GenParticle(ap,APartID,1);
HepMC::GenParticle* APart = new HepMC::GenParticle(ap,APartID,1);
APart->suggest_barcode( barcode ) ;
barcode++ ;
Vtx->add_particle_out(APart) ;
Expand Down
13 changes: 13 additions & 0 deletions RecoMuon/CosmicMuonProducer/interface/HeaderForQC8.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,13 @@
#ifndef CosmicMuonProducer_HeaderForQC8_H
#define CosmicMuonProducer_HeaderForQC8_H


#define QC8FLAG_SEEDINFO_REFVERTROLL18 0x00000001
#define QC8FLAG_SEEDINFO_DIFFCOL1 0x00000002
#define QC8FLAG_SEEDINFO_DIFFCOL2 0x00000004


#endif //CosmicMuonProducer_HeaderForQC8_H



2 changes: 1 addition & 1 deletion RecoMuon/CosmicMuonProducer/src/GEMCosmicMuon.cc
Original file line number Diff line number Diff line change
Expand Up @@ -292,7 +292,7 @@ Trajectory GEMCosmicMuon::makeTrajectory(TrajectorySeed seed, MuonTransientTrack
//const DetLayer* layer = theService->detLayerGeometry()->idToLayer( ch->id().rawId() );
std::shared_ptr<MuonTransientTrackingRecHit> tmpRecHit;
tsosCurrent = theService->propagator("SteppingHelixPropagatorAny")->propagate(tsosCurrent,ch->surface());
if (!tsosCurrent.isValid()) continue;
if (!tsosCurrent.isValid()) return Trajectory();
GlobalPoint tsosGP = tsosCurrent.freeTrajectoryState()->position();
//TrackingRecHit *tmpRecHit = new TrackingRecHit(ch);
//cout << "tsosGP "<< tsosGP <<endl;
Expand Down
Loading