Go to the documentation of this file.00001
00002
00003
00004
00005
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022 #include <memory>
00023 #include <iostream>
00024 #include <fstream>
00025 #include <vector>
00026 #include <string>
00027
00028 #include "FWCore/Framework/interface/Frameworkfwd.h"
00029 #include "FWCore/Framework/interface/EDAnalyzer.h"
00030
00031 #include "FWCore/Framework/interface/Event.h"
00032 #include "FWCore/Framework/interface/EventSetup.h"
00033 #include "FWCore/Framework/interface/MakerMacros.h"
00034 #include "FWCore/Framework/interface/ESHandle.h"
00035 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00036
00037 #include "CommonTools/UtilAlgos/interface/TFileService.h"
00038 #include "FWCore/ServiceRegistry/interface/Service.h"
00039
00040 #include "DataFormats/Candidate/interface/Candidate.h"
00041 #include "DataFormats/RecoCandidate/interface/RecoCaloTowerCandidate.h"
00042 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
00043 #include "DataFormats/Common/interface/EDProduct.h"
00044 #include "DataFormats/Common/interface/Ref.h"
00045
00046 #include "CondFormats/DataRecord/interface/HeavyIonRPRcd.h"
00047 #include "CondFormats/HIObjects/interface/CentralityTable.h"
00048 #include "CondCore/DBOutputService/interface/PoolDBOutputService.h"
00049 #include "RecoHI/HiEvtPlaneAlgos/interface/HiEvtPlaneList.h"
00050
00051 #include "CondFormats/HIObjects/interface/RPFlatParams.h"
00052 #include "TFile.h"
00053 #include "TH1.h"
00054 #include "TH2D.h"
00055 #include "TH2F.h"
00056 #include "TTree.h"
00057 #include "TH1I.h"
00058 #include "TF1.h"
00059 #include "TList.h"
00060 #include "TString.h"
00061 #include <time.h>
00062 #include <cstdlib>
00063 #include <iostream>
00064 #include <vector>
00065 using namespace std;
00066 using namespace hi;
00067
00068
00069
00070
00071
00072 class MoveFlatParamsToDB : public edm::EDAnalyzer {
00073 public:
00074 explicit MoveFlatParamsToDB(const edm::ParameterSet&);
00075 ~MoveFlatParamsToDB();
00076
00077
00078 private:
00079 virtual void beginJob() ;
00080 virtual void analyze(const edm::Event&, const edm::EventSetup&);
00081 virtual void endJob() ;
00082
00083 TFile * inFile;
00084
00085 TH1D * x[NumEPNames];
00086 TH1D * y[NumEPNames];
00087 TH1D * xycnt[NumEPNames];
00088 string rpname[NumEPNames];
00089 int RPNameIndx[NumEPNames];
00090 int RPSubEvnt[NumEPNames];
00091 RPFlatParams * rpFlat;
00092 int nRP;
00093
00094
00095 };
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108 MoveFlatParamsToDB::MoveFlatParamsToDB(const edm::ParameterSet& iConfig)
00109
00110 {
00111 cout<<"Enter MoveFlatParamsToDB"<<endl;
00112
00113 inFile = new TFile("data/rpflat_combined.root");
00114 if(inFile->IsZombie()) {
00115 cout<<"file not found"<<endl;
00116 }
00117 TList * list = ((TDirectory *)inFile->Get("hiEvtPlaneFlatCalib"))->GetListOfKeys();
00118 int indx =0;
00119 int cnt = 0;
00120 for(int i = 0; i<NumEPNames; i++) {
00121 x[i]=0;
00122 y[i]=0;
00123 xycnt[i]=0;
00124 }
00125 while(indx >=0 && indx<NumEPNames) {
00126 int EPNameIndx = -1;
00127 TString name = list->At(indx)->GetName();
00128 if(!name.Contains("cent")&&!name.Contains("vtx")&&!name.Contains("MidEtaTrackRescor")) {
00129 for(int i = 0; i<NumEPNames; i++) {
00130 if(name.CompareTo(EPNames[i])==0) {
00131 EPNameIndx = i;
00132 break;
00133 }
00134 }
00135 if(EPNameIndx <0) cout<<"A bad reaction plane name has been encountered: "<<name.Data()<<endl;
00136 RPNameIndx[cnt]=EPNameIndx;
00137 TString sname = name;
00138 x[cnt] = (TH1D *) inFile->Get(Form("hiEvtPlaneFlatCalib/%s/x_%s",name.Data(),name.Data()));
00139 y[cnt] = (TH1D *) inFile->Get(Form("hiEvtPlaneFlatCalib/%s/y_%s",name.Data(),name.Data()));
00140 xycnt[cnt] = (TH1D *) inFile->Get(Form("hiEvtPlaneFlatCalib/%s/cnt_%s",name.Data(),name.Data()));
00141 rpname[cnt]=sname;
00142 if(!x[cnt]) cout<<"bad x"<<endl;
00143 if(!y[cnt]) cout<<"bad y"<<endl;
00144 if(!xycnt[cnt]) cout<<"bad cnt"<<endl;
00145 if(x[cnt] && xycnt[cnt] && y[cnt]) {
00146 x[cnt]->Divide(xycnt[cnt]);
00147 y[cnt]->Divide(xycnt[cnt]);
00148 }
00149 ++cnt;
00150 if(cnt>NumEPNames||cnt>50) {
00151 cout<<"Maximum number of reaction planes exceeded!"<<endl;
00152 break;
00153 }
00154 }
00155
00156 if(list->At(indx)==list->Last())
00157 indx = -1;
00158 else
00159 ++indx;
00160 }
00161 nRP = cnt;
00162 cout<<"nRP = "<<nRP<<endl;
00163 }
00164
00165
00166 MoveFlatParamsToDB::~MoveFlatParamsToDB()
00167 {
00168
00169
00170
00171
00172 }
00173
00174
00175
00176
00177
00178
00179
00180 void
00181 MoveFlatParamsToDB::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00182 {
00183 using namespace edm;
00184 rpFlat = new RPFlatParams();
00185 rpFlat->m_table.reserve(x[0]->GetNbinsX());
00186 cout<<"Size of table: "<<x[0]->GetNbinsX()<<endl;
00187 for(int j = 0; j<x[0]->GetNbinsX();j++) {
00188 RPFlatParams::EP * thisBin = new RPFlatParams::EP();
00189 for(int i = 0; i<nRP; i++) {
00190 thisBin->x[i] = x[i]->GetBinContent(j+1);
00191 thisBin->y[i] = y[i]->GetBinContent(j+1);
00192 thisBin->RPNameIndx[i]=RPNameIndx[i];
00193 }
00194 rpFlat->m_table.push_back(*thisBin);
00195 if(thisBin) delete thisBin;
00196 }
00197 cout<<"Number of RP: "<<nRP<<endl;
00198 edm::Service<cond::service::PoolDBOutputService> poolDbService;
00199 if(poolDbService.isAvailable())
00200 poolDbService->writeOne( rpFlat,poolDbService->beginOfTime(),"HeavyIonRPRcd");
00201 cout<<"DONE"<<endl;
00202 }
00203
00204
00205
00206 void
00207 MoveFlatParamsToDB::beginJob()
00208 {
00209 }
00210
00211
00212 void
00213 MoveFlatParamsToDB::endJob() {
00214 }
00215
00216
00217 DEFINE_FWK_MODULE(MoveFlatParamsToDB);