Fixing error in Geometry, and corresponding codes to plot on TH2Poly

Geometry/src/HGCalTBCellVertices.cc

@@ -56,7 +56,7 @@ std::pair<double, double> HGCalTBCellVertices::GetCellCentreCoordinates(int laye
 	bool ValidFlag   = Top.iu_iv_valid(layer, sensor_iu, sensor_iv, iu, iv, sensorsize);
 	if(ValidFlag) {    
 		centre_x_tmp = iu * x0 + iv * vx0 + sensor_iu*X0 + sensor_iv*VX0;
-		centre_y_tmp = iv * vy0 + iv * vy0 + sensor_iv*VY0;
+		centre_y_tmp = iv * vy0 + sensor_iv*VY0;
 		auto point = RotateLayer(std::make_pair(centre_x_tmp, centre_y_tmp), TEST_BEAM_LAYER_ROTATION);
 //		if(flipX==true) point.first = - point.first;
 		return point;

Geometry/src/HGCalTBTopology.cc

@@ -11,13 +11,14 @@ bool HGCalTBTopology::iu_iv_valid(int layer, int sensor_iu, int sensor_iv, int i
         int aiv = abs(iv);
         int iuc = (iv < 0) ? (-iu) : (iu);
 
-        int aIv = abs(sensor_iv);
-        int Iuc = (sensor_iv < 0)?(-sensor_iu):(sensor_iu);
+	int aIv = sensor_iv;
+        int Iuc = sensor_iu;
 
-        if(layer <= 28) Is_Valid_sensor_iu_iv = (Iuc == 0 && aIv == 0);
+	if(layer <= 4) Is_Valid_sensor_iu_iv = (Iuc == 0 && aIv == 0);
         else{
-                if(aIv == 0) Is_Valid_sensor_iu_iv = (Iuc >= -1 && Iuc <= 1);
-                else if(aIv == 1) Is_Valid_sensor_iu_iv = (Iuc >= -1 && Iuc <= 0);
+                if((Iuc == 1) && (aIv == -1) ) Is_Valid_sensor_iu_iv = 1;
+                else if((Iuc == 0) && (aIv == 0) ) Is_Valid_sensor_iu_iv = 1;
+                else if  ( (Iuc == 0) && (aIv == -1) ) Is_Valid_sensor_iu_iv = 1;
         }
 
 	if(layer <= 40 && Is_Valid_sensor_iu_iv) {

RawToDigi/plugins/HGCalTBRawDataSource.cc

@@ -9,6 +9,8 @@
 #include <ctime>
 #include <cmath>
 
+using namespace std;
+
 HGCalTBRawDataSource::HGCalTBRawDataSource(const edm::ParameterSet & pset, edm::InputSourceDescription const& desc) :
   edm::ProducerSourceFromFiles(pset, desc, true),
   m_electronicMap(pset.getUntrackedParameter<std::string>("ElectronicMap","HGCal/CondObjects/data/map_CERN_Hexaboard_28Layers.txt")),
@@ -202,6 +204,7 @@ void HGCalTBRawDataSource::produce(edm::Event & e)
       else{ 
 	HGCalTBDetId did = essource_.emap_.eid2detId(eid);
 	detids.push_back(did);
+//	cout<<endl<<" ski = "<<iski<<"  Chan = "<<ichan<<" Layer = "<<did.layer()<<" Sensor IU = "<< did.sensorIU()<<" Sensor IV = "<< did.sensorIV()<<" iu = "<< did.iu()<<" iv = "<<did.iv()<<endl;
       }
     }
     std::vector<uint16_t> vdata;vdata.insert(vdata.end(),m_decodedData.at(iski).begin(),m_decodedData.at(iski).end());

RawToDigi/plugins/PedestalPlotter.cc

@@ -23,6 +23,8 @@
 #include <iomanip>
 #include <set>
 
+using namespace std;
+
 struct hgcal_channel{
   hgcal_channel() : key(0),
 		    medianHG(0.),
@@ -103,10 +105,15 @@ void PedestalPlotter::analyze(const edm::Event& event, const edm::EventSetup& se
   edm::Handle<HGCalTBSkiroc2CMSCollection> skirocs;
   event.getByToken(m_HGCalTBSkiroc2CMSCollection, skirocs);
 
+
   m_numberOfBoards = skirocs->size()/HGCAL_TB_GEOMETRY::N_SKIROC_PER_HEXA;
 
   m_evtID++;
 
+
+//  if(m_evtID > 100) return;
+
+
   if( !skirocs->size() ) return;
 
   for( size_t iski=0;iski<skirocs->size(); iski++ ){
@@ -121,6 +128,9 @@ void PedestalPlotter::analyze(const edm::Event& event, const edm::EventSetup& se
 	setOfConnectedDetId.insert(p);
       }
       else continue;
+
+//      cout<<endl<<" Ski = "<<iski<<"  Chan = "<<ichan<<" Layer = "<<detid.layer()<<" Sensor IU = "<< detid.sensorIU()<<" Sensor IV = "<< detid.sensorIV()<<" iu = "<< detid.iu()<<" iv = "<<detid.iv()<<endl;
+
       for( size_t it=0; it<NUMBER_OF_SCA; it++ ){
 	if( rollpositions[it]<2 ){ //keep only 2 first time sample for pedestal evaluation
 	  uint32_t key=iboard*100000+(iski%HGCAL_TB_GEOMETRY::N_SKIROC_PER_HEXA)*10000+ichan*100+it;
@@ -182,16 +192,32 @@ void PedestalPlotter::endJob()
   std::map<int,TH2Poly*>  chanMap;
   std::ostringstream os( std::ostringstream::ate );
   TH2Poly *h;
+  int Board_IU = 0;
+  int Board_IV = 0;	
+  int Board_Layer = 0;
   for(size_t ib = 0; ib<m_numberOfBoards; ib++) {
+    if( (ib == 5) || (ib == 8) ){
+	Board_IU = 0;
+	Board_IV = -1;
+    }
+    if( (ib == 4) || (ib == 7) ){
+        Board_IU = 1;
+        Board_IV = -1;
+    }
+
+    if(ib <= 3) Board_Layer = ib + 1;
+    else if( (ib == 4) || (ib == 5) || (ib == 6) ) Board_Layer = 5;
+    else if( (ib == 7) || (ib == 8) || (ib == 9) ) Board_Layer = 6;
+
     os.str("");
     os << "HexaBoard" << ib ;
     TFileDirectory dir = fs->mkdir( os.str().c_str() );
     h=dir.make<TH2Poly>();
     os.str("");
-    os<<"ChannelMapping";
+    os<<"ChannelMapping_"<<"Layer_"<<Board_Layer<<"_Sensor_IU_"<<Board_IU<<"_Sensor_IV_"<<Board_IV;
     h->SetName(os.str().c_str());
     h->SetTitle(os.str().c_str());
-    InitTH2Poly(*h, (int)ib, 0, 0);
+    InitTH2Poly(*h, Board_Layer, Board_IU, Board_IV);
     chanMap.insert( std::pair<int,TH2Poly*>(ib,h) );
     TFileDirectory hgpdir = dir.mkdir( "HighGainPedestal" );
     TFileDirectory lgpdir = dir.mkdir( "LowGainPedestal" );
@@ -201,34 +227,35 @@ void PedestalPlotter::endJob()
 
       h=hgpdir.make<TH2Poly>();
       os.str("");
-      os<<"SCA"<<it;
+      os<<"Layer_"<<Board_Layer<<"_Sensor_IU_"<<Board_IU<<"_Sensor_IV_"<<Board_IV<<"_SCA_"<<it;
       h->SetName(os.str().c_str());
       h->SetTitle(os.str().c_str());
-      InitTH2Poly(*h, (int)ib, 0, 0);
+      InitTH2Poly(*h, Board_Layer, Board_IU, Board_IV);
       hgMeanMap.insert( std::pair<int,TH2Poly*>(100*ib+it,h) );
 
       h=lgpdir.make<TH2Poly>();
       os.str("");
+
       os<<"SCA"<<it;
       h->SetName(os.str().c_str());
       h->SetTitle(os.str().c_str());
-      InitTH2Poly(*h, (int)ib, 0, 0);
+      InitTH2Poly(*h, Board_Layer, Board_IU, Board_IV);
       lgMeanMap.insert( std::pair<int,TH2Poly*>(100*ib+it,h) );
 
       h=hgndir.make<TH2Poly>();
       os.str("");
       os<<"SCA"<<it;
       h->SetName(os.str().c_str());
       h->SetTitle(os.str().c_str());
-      InitTH2Poly(*h, (int)ib, 0, 0);
+      InitTH2Poly(*h, Board_Layer, Board_IU, Board_IV);
       hgRMSMap.insert( std::pair<int,TH2Poly*>(100*ib+it,h) );
 
       h=lgndir.make<TH2Poly>();
       os.str("");
       os<<"SCA"<<it;
       h->SetName(os.str().c_str());
       h->SetTitle(os.str().c_str());
-      InitTH2Poly(*h, (int)ib, 0, 0);
+      InitTH2Poly(*h, Board_Layer, Board_IU, Board_IV);
       lgRMSMap.insert( std::pair<int,TH2Poly*>(100*ib+it,h) );
 
     }
@@ -266,25 +293,29 @@ void PedestalPlotter::endJob()
     }
     HGCalTBDetId detid=(*it).second;
     CellCentreXY = TheCell.GetCellCentreCoordinatesForPlots( detid.layer(), detid.sensorIU(), detid.sensorIV(), detid.iu(), detid.iv(), m_sensorsize );
+
     double iux = (CellCentreXY.first < 0 ) ? (CellCentreXY.first + HGCAL_TB_GEOMETRY::DELTA) : (CellCentreXY.first - HGCAL_TB_GEOMETRY::DELTA) ;
     double iuy = (CellCentreXY.second < 0 ) ? (CellCentreXY.second + HGCAL_TB_GEOMETRY::DELTA) : (CellCentreXY.second - HGCAL_TB_GEOMETRY::DELTA);
+
+//    cout<<endl<<" iboard = "<<iboard<<" IU = "<<detid.sensorIU()<<" IV = "<<detid.sensorIV()<<" X = "<<iux<<" Y = "<<iuy<<endl;
+
     for( size_t it=0; it<NUMBER_OF_SCA; it++ ){
       int key=iboard*100000+iski*10000+ichan*100+it;
       std::map<int,hgcal_channel>::iterator iter=m_channelMap.find(key);
       float hgMean=iter->second.medianHG;
       float lgMean=iter->second.medianLG;
       float hgRMS=iter->second.rmsHG;
       float lgRMS=iter->second.rmsLG;
-      hgMeanMap[ 100*iboard+it ]->Fill(iux/2 , iuy, hgMean );
-      lgMeanMap[ 100*iboard+it ]->Fill(iux/2 , iuy, lgMean );
-      hgRMSMap[ 100*iboard+it ]->Fill(iux/2 , iuy, hgRMS );
-      lgRMSMap[ 100*iboard+it ]->Fill(iux/2 , iuy, lgRMS );
+      hgMeanMap[ 100*iboard+it ]->Fill(iux , iuy, hgMean );
+      lgMeanMap[ 100*iboard+it ]->Fill(iux , iuy, lgMean );
+      hgRMSMap[ 100*iboard+it ]->Fill(iux , iuy, hgRMS );
+      lgRMSMap[ 100*iboard+it ]->Fill(iux , iuy, lgRMS );
       if( m_writePedestalFile ){
 	pedestalHG << " " << hgMean << " " << hgRMS;
 	pedestalLG << " " << lgMean << " " << lgRMS;;
       }
     }
-    chanMap[ iboard ]->Fill(iux/2 , iuy, iski*1000+ichan );
+    chanMap[ iboard ]->Fill(iux , iuy, iski*1000+ichan );
     if( m_writePedestalFile ){
       pedestalHG << std::endl;
       pedestalLG << std::endl;

RawToDigi/plugins/RawDataPlotter.cc

@@ -168,7 +168,7 @@ void RawDataPlotter::analyze(const edm::Event& event, const edm::EventSetup& set
 	  CellCentreXY = TheCell.GetCellCentreCoordinatesForPlots( detid.layer(), detid.sensorIU(), detid.sensorIV(), detid.iu(), detid.iv(), m_sensorsize );
 	  double iux = (CellCentreXY.first < 0 ) ? (CellCentreXY.first + HGCAL_TB_GEOMETRY::DELTA) : (CellCentreXY.first - HGCAL_TB_GEOMETRY::DELTA) ;
 	  double iuy = (CellCentreXY.second < 0 ) ? (CellCentreXY.second + HGCAL_TB_GEOMETRY::DELTA) : (CellCentreXY.second - HGCAL_TB_GEOMETRY::DELTA);
-	  polyMap[ 100*iboard+it ]->Fill(iux/2 , iuy, skiroc.ADCHigh(ichan,it) );
+	  polyMap[ 100*iboard+it ]->Fill(iux , iuy, skiroc.ADCHigh(ichan,it) );
 	}
       }
     }

RawToDigi/plugins/RawHitPlotter.cc

@@ -254,7 +254,7 @@ void RawHitPlotter::analyze(const edm::Event& event, const edm::EventSetup& setu
       CellCentreXY=TheCell.GetCellCentreCoordinatesForPlots(hit.detid().layer(),hit.detid().sensorIU(),hit.detid().sensorIV(),hit.detid().iu(),hit.detid().iv(),m_sensorsize);
       double iux = (CellCentreXY.first < 0 ) ? (CellCentreXY.first + delta) : (CellCentreXY.first - delta) ;
       double iuy = (CellCentreXY.second < 0 ) ? (CellCentreXY.second + delta) : (CellCentreXY.second - delta);
-      polyMap[ 100*iboard+it ]->Fill(iux/2 , iuy, highGain);
+      polyMap[ 100*iboard+it ]->Fill(iux , iuy, highGain);
     }
   }
 }
@@ -347,10 +347,10 @@ void RawHitPlotter::endJob()
       float lgMean=cif->meanLGMap[it]/cif->counterLGMap[it];
       float hgRMS=std::sqrt(cif->rmsHGMap[it]/cif->counterHGMap[it]-cif->meanHGMap[it]/cif->counterHGMap[it]*cif->meanHGMap[it]/cif->counterHGMap[it]);
       float lgRMS=std::sqrt(cif->rmsLGMap[it]/cif->counterLGMap[it]-cif->meanLGMap[it]/cif->counterLGMap[it]*cif->meanLGMap[it]/cif->counterLGMap[it]);
-      hgMeanMap[ 100*iboard+it ]->Fill(iux/2 , iuy, hgMean );
-      lgMeanMap[ 100*iboard+it ]->Fill(iux/2 , iuy, lgMean );
-      hgRMSMap[ 100*iboard+it ]->Fill(iux/2 , iuy, hgRMS );
-      lgRMSMap[ 100*iboard+it ]->Fill(iux/2 , iuy, lgRMS );
+      hgMeanMap[ 100*iboard+it ]->Fill(iux , iuy, hgMean );
+      lgMeanMap[ 100*iboard+it ]->Fill(iux , iuy, lgMean );
+      hgRMSMap[ 100*iboard+it ]->Fill(iux , iuy, hgRMS );
+      lgRMSMap[ 100*iboard+it ]->Fill(iux , iuy, lgRMS );
     }
   }
 }

Reco/plugins/HGCalTBRecHitProducer.cc

@@ -56,7 +56,7 @@ void HGCalTBRecHitProducer::produce(edm::Event& event, const edm::EventSetup& iS
 
     std::vector<double> sampleHG, sampleLG, sampleT;
 
-    float highGain, lowGain, totGain;
+    float highGain = -500, lowGain = -500, totGain;
     int hgStatus = -1;
     int lgStatus = -1;
     float timeHG = 0.;

Reco/plugins/RecHitPlotter.cc

@@ -187,8 +187,8 @@ void RecHitPlotter::analyze(const edm::Event& event, const edm::EventSetup& setu
       CellCentreXY=TheCell.GetCellCentreCoordinatesForPlots(hit.id().layer(),hit.id().sensorIU(),hit.id().sensorIV(),hit.id().iu(),hit.id().iv(),m_sensorsize);
       double iux = (CellCentreXY.first < 0 ) ? (CellCentreXY.first + HGCAL_TB_GEOMETRY::DELTA) : (CellCentreXY.first - HGCAL_TB_GEOMETRY::DELTA) ;
       double iuy = (CellCentreXY.second < 0 ) ? (CellCentreXY.second + HGCAL_TB_GEOMETRY::DELTA) : (CellCentreXY.second - HGCAL_TB_GEOMETRY::DELTA);
-      polyMapHG[ (eid.iskiroc()-1)/HGCAL_TB_GEOMETRY::N_SKIROC_PER_HEXA ]->Fill(iux/2 , iuy, hit.energyHigh());
-      polyMapLG[ (eid.iskiroc()-1)/HGCAL_TB_GEOMETRY::N_SKIROC_PER_HEXA ]->Fill(iux/2 , iuy, hit.energyLow());
+      polyMapHG[ (eid.iskiroc()-1)/HGCAL_TB_GEOMETRY::N_SKIROC_PER_HEXA ]->Fill(iux , iuy, hit.energyHigh());
+      polyMapLG[ (eid.iskiroc()-1)/HGCAL_TB_GEOMETRY::N_SKIROC_PER_HEXA ]->Fill(iux , iuy, hit.energyLow());
     }
   }
   m_h_hgSum->Fill( energyHighSum );

runPedestalEval.sh

@@ -0,0 +1,43 @@
+#!/bin/bash
+
+run=$1
+
+path="/afs/cern.ch/work/r/rchatter/TestBeam_July/CMSSW_8_0_21/src/HGCal"
+CFG="${path}/unpack2017_cfg.py"
+TOP=${PWD}
+
+echo ${path}
+echo ${CFG}
+
+cp ${CFG} ./cfg.py
+
+cd ${path}
+eval `scramv1 runtime -sh` 
+cd ${TOP}
+
+eosPath="/eos/cms/store/group/dpg_hgcal/tb_hgcal/july2017/July2017_TB_data_orm"
+echo "Download ${eosPath}/HexaData_Run${run}.raw"
+
+# THIS IS DO DOWNLOAD THE TIMING TXT FILES WHICH YOU DON'T NEED
+ echo "Download ${eosPath}/HexaData_Run${run}_TIMING_RDOUT_ORM0.txt"
+ xrdcp -f root://eoscms.cern.ch/${eosPath}/HexaData_Run${run}_TIMING_RDOUT_ORM0.txt HexaData_Run${run}_TIMING_RDOUT_ORM0.txt
+ echo "Download ${eosPath}/HexaData_Run${run}_TIMING_RDOUT_ORM1.txt"
+ xrdcp -f root://eoscms.cern.ch/${eosPath}/HexaData_Run${run}_TIMING_RDOUT_ORM1.txt HexaData_Run${run}_TIMING_RDOUT_ORM1.txt
+ echo "Download ${eosPath}/HexaData_Run${run}_TIMING_RDOUT_ORM2.txt"
+ xrdcp -f root://eoscms.cern.ch/${eosPath}/HexaData_Run${run}_TIMING_RDOUT_ORM2.txt HexaData_Run${run}_TIMING_RDOUT_ORM2.txt
+
+cmsRun unpack2017_cfg.py runNumber=${run} dataFolder="./" outputFolder="./"
+
+# THIS IS TO SAVE OUTPUT ON EOS : PLEASE DON'T DO THAT
+# eosOutputPath="/eos/cms/store/group/dpg_hgcal/tb_hgcal/july2017/pedestalFiles"
+# 
+# echo "Upload ${eosOutputPath}/PedestalOutput_${run}.root"
+# xrdcp -f PedestalOutput_${run}.root root://eoscms.cern.ch/${eosOutputPath}/PedestalOutput_${run}.root
+# echo "Upload ${eosOutputPath}/pedestalHG_${run}.txt"
+# xrdcp -f pedestalHG_${run}.txt root://eoscms.cern.ch/${eosOutputPath}/pedestalHG_${run}.txt
+# echo "Upload ${eosOutputPath}/pedestalLG_${run}.txt"
+# xrdcp -f pedestalLG_${run}.txt root://eoscms.cern.ch/${eosOutputPath}/pedestalLG_${run}.txt
+# 
+# eosOutputPath="/eos/cms/store/group/dpg_hgcal/tb_hgcal/july2017/HGCalTBSkiroc2CMS"
+# echo "Upload ${eosOutputPath}/cmsswEvents_Run${run}.root"
+# xrdcp -f cmsswEvents.root root://eoscms.cern.ch/${eosOutputPath}/cmsswEvents_Run${run}.root

unpack2017_cfg.py

@@ -31,7 +31,7 @@
 if not os.path.isdir(options.dataFolder):
     sys.exit("Error: Data folder not found or inaccessible!")
 
-electronicMap="HGCal/CondObjects/data/map_CERN_Hexaboard_28Layers_AllFlipped.txt"
+electronicMap="HGCal/CondObjects/data/map_CERN_Hexaboard_July_6Layers.txt"
 pedestalHighGain=options.outputFolder+"/pedestalHG_"+str(options.runNumber)+".txt"
 pedestalLowGain=options.outputFolder+"/pedestalLG_"+str(options.runNumber)+".txt"
 
@@ -53,7 +53,7 @@
                             NumberOfBytesForTheTrailer=cms.untracked.uint32(4),
                             NumberOfBytesForTheEventTrailers=cms.untracked.uint32(4),
                             NSkipEvents=cms.untracked.uint32(0),
-                            ReadTXTForTiming=cms.untracked.bool(True),
+                            ReadTXTForTiming=cms.untracked.bool(False),
                             timingFiles=cms.vstring("%s/HexaData_Run%04d_TIMING_RDOUT_ORM0.txt"%(options.dataFolder,options.runNumber),
                                                     "%s/HexaData_Run%04d_TIMING_RDOUT_ORM1.txt"%(options.dataFolder,options.runNumber),
                                                     "%s/HexaData_Run%04d_TIMING_RDOUT_ORM2.txt"%(options.dataFolder,options.runNumber))