Warning, /Validation/DTRecHits/test/writeSummaryTable.r is written in an unsupported language. File is not indexed.
0001 //------------------------------
0002 //
0003 // Write table of true hit resolutions, pull width, and efficiencies from DQM files produced in local mode
0004 //
0005 // Usage:
0006 // .x writeSummaryTable.r("inputFile.root",doEff,doRes,doPull)
0007 //
0008 // Author: N. Amapane
0009 //
0010 //------------------------------
0011 #include <iomanip>
0012
0013
0014 void writeSummaryTable(){
0015 cout << endl << "Usage: .x writeSummaryTable.r(\"inputFile.root\",doEff,doEffab,doRes,doResab,doPull)" << endl << endl;
0016 }
0017
0018
0019 void writeSummaryTable(TString filename, bool doEff=true, bool doEffab = true, bool doRes=true, bool doResab = true, bool doPull=true,bool doPullabxy=true, bool doSegRes=true, bool doNSeg= true, bool doPar = false) {
0020
0021 //----------------------------------------------------------------------
0022 // Configurable options
0023
0024 // Interval for the fit of residual distributions
0025 float nsigma = 2.;
0026
0027 //----------------------------------------------------------------------
0028
0029 if (! TString(gSystem->GetLibraries()).Contains("Histograms_h")) {
0030 gROOT->LoadMacro("$CMSSW_BASE/src/Validation/DTRecHits/interface/Histograms.h+");
0031 gROOT->LoadMacro("macros.C");
0032 gROOT->LoadMacro("ranges.C");
0033 }
0034
0035 float cmToMicron = 10000.;
0036
0037 TFile *file = new TFile(filename);
0038
0039 TString table = filename.Replace(filename.Length()-5,5,"_summary.txt");
0040 ofstream f(table,ios_base::out);
0041 f << fixed;
0042 f << "# W St sec SL effS1RPhi effS3RPhi effSeg resHit pullHit meanAngle sigmaAngle meanSegPos sigmaSegPos" << endl;
0043
0044 // All sectors are collapsed together as of now
0045 int smin = 0;
0046 int smax = 0;
0047
0048 for (int wheel = -2; wheel<3; ++wheel) {
0049 for (int station =1; station<=4; ++station) {
0050 for (int sector = smin; sector<=smax; ++sector) {
0051
0052 if (station!=4 && sector>12) continue;
0053
0054 double stheta = -1.;
0055 double sphi = -1.;
0056 double ptheta = -1.;
0057 double pphi = -1.;
0058 double salpha_res = -1.;
0059 double sbeta_res = -1.;
0060 double malpha_res = -1.;
0061 double mbeta_res = -1.;
0062 double salpha_pull = -1.;
0063 double sbeta_pull = -1.;
0064 double malpha_pull = -1.;
0065 double mbeta_pull = -1.;
0066 double sx = -1.;
0067 double sy = -1.;
0068 double mx = -1.;
0069 double my = -1.;
0070 double sx_pull = -1.;
0071 double sy_pull = -1.;
0072 double mx_pull = -1.;
0073 double my_pull = -1.;
0074
0075 float effS1RPhi=0.;
0076 float effS3RPhi=0.;
0077 float effS1RZ=0.;
0078 float effS3RZ=0.;
0079 float effSeg = 0;
0080
0081
0082
0083 HRes1DHit *hResPhi1 = new HRes1DHit(file, wheel, station, 1, "S3");
0084 HRes1DHit *hResTheta = new HRes1DHit(file, wheel, station, 2, "S3");
0085
0086 HEff1DHit* hEffS1RPhi= new HEff1DHit(file, wheel, station, 1, "S1");
0087 HEff1DHit* hEffS3RPhi= new HEff1DHit(file, wheel, station, 1, "S3");
0088 HEff1DHit* hEffS1RZ=0;
0089 HEff1DHit* hEffS3RZ=0;
0090
0091 HRes4DHit* hRes4D= new HRes4DHit(file, wheel, station, 0);
0092 HEff4DHit* hEff4D = new HEff4DHit(file, wheel, station, 0);
0093
0094
0095 if (station!=4) {
0096 hEffS1RZ= new HEff1DHit(file, wheel, station, 2, "S1");
0097 hEffS3RZ= new HEff1DHit(file, wheel, station, 2, "S3");
0098 }
0099
0100
0101
0102 TCanvas c;
0103
0104 //--- Hit efficiency
0105 if (doEff) {
0106 if (station!=4) {
0107 effS1RZ = hEffS1RZ->hDistRecHit->Integral()/hEffS1RZ->hDistMuSimHit->Integral();
0108 effS3RZ = hEffS3RZ->hDistRecHit->Integral()/hEffS3RZ->hDistMuSimHit->Integral();
0109 }
0110 effS1RPhi = hEffS1RPhi->hDistRecHit->Integral()/hEffS1RPhi->hDistMuSimHit->Integral();
0111 effS3RPhi = hEffS3RPhi->hDistRecHit->Integral()/hEffS3RPhi->hDistMuSimHit->Integral();
0112
0113 cout << " " << wheel << " " << station << " " << sector << " effPhi: " << effS1RPhi << " " << effS3RPhi;
0114 if (station!=4) cout << " effZ: " << effS1RZ << " " << effS3RZ;
0115 cout << endl;
0116 }
0117
0118
0119 //--- Segment efficiency
0120 if (doEffab) {
0121 if (station!=4) {
0122
0123 effSeg = hEff4D->hBetaRecHit->Integral()/hEff4D->hBetaSimSegm->Integral();
0124
0125 cout << " " << wheel << " " << station << " " << sector << " effSeg: " << effSeg << endl;
0126
0127
0128 }
0129 }
0130
0131
0132 //--- Hit resolution
0133 if (doRes) {
0134 if (station!=4) {
0135 TH1F* tmpTheta = (TH1F*) hResTheta->hRes->Clone("tmpTheta");
0136 tmpTheta->Rebin(2);
0137 TF1* ftheta = drawGFit(tmpTheta, nsigma, -2. , 2.);
0138 stheta = ftheta->GetParameter("Sigma");
0139 }
0140 TH1F* tmpPhi = (TH1F*) hResPhi1->hRes->Clone("tmpPhi");
0141 tmpPhi->Rebin(2);
0142 TF1* fphi = drawGFit(tmpPhi, nsigma, -2. , 2. );
0143 sphi = fphi->GetParameter("Sigma");
0144
0145 cout << " " << wheel << " " << station << " " << sector << " sphi: " << sphi;
0146 if (station!=4) cout << " stheta: " << stheta ;
0147 cout << endl;
0148 }
0149
0150
0151 //--- Hit pull
0152 if (doPull) {
0153 if (station!=4) {
0154 TH1F* tmpTheta = (TH1F*) hResTheta->hPull->Clone("tmpTheta");
0155 TF1* ftheta = drawGFit(tmpTheta, nsigma, -2. , 2.);
0156 ptheta = ftheta->GetParameter("Sigma");
0157 }
0158 TH1F* tmpPhi = (TH1F*) hResPhi1->hPull->Clone("tmpPhi");
0159 TF1* fphi = drawGFit(tmpPhi, nsigma, -2. , 2. );
0160 pphi = fphi->GetParameter("Sigma");
0161
0162 cout << " " << wheel << " " << station << " " << sector
0163 << " pphi: " << pphi;
0164 if (station!=4) {
0165 cout << " ptheta: " << ptheta ;
0166 }
0167 cout << endl;
0168 }
0169
0170
0171 //--- Segment position resolution
0172 if (doSegRes) {
0173 if (station!=4) {
0174 TH1F* tmpY = (TH1F*) hRes4D->hResYRZ->Clone("tmpY");
0175 TF1* fy = drawGFit(tmpY, nsigma, -2. , 2. );
0176 sy = fy->GetParameter("Sigma");
0177 my = fy->GetParameter("Mean");
0178 }
0179
0180 TH1F* tmpX = (TH1F*) hRes4D->hResX->Clone("tmpX");
0181 TF1* fx = drawGFit(tmpX, nsigma, -2. , 2.);
0182 sx = fx->GetParameter("Sigma");
0183 mx = fx->GetParameter("Mean");
0184 }
0185
0186 //--- Segment angle resolution
0187 if (doResab) {
0188 if (station!=4) {
0189 TH1F* tmpBeta = (TH1F*) hRes4D->hResBeta->Clone("tmpBeta");
0190 tmpBeta->Rebin(2);
0191 TF1* fbeta = drawGFit(tmpBeta, nsigma, -2. , 2.);
0192 sbeta_res = fbeta->GetParameter("Sigma");
0193 mbeta_res = fbeta->GetParameter("Mean");
0194 }
0195
0196 TH1F* tmpAlpha = (TH1F*) hRes4D->hResAlpha->Clone("tmpAlpha");
0197 // tmpAlpha->Rebin(2);
0198 TF1* falpha = drawGFit(tmpAlpha, nsigma, -2. , 2. );
0199 salpha_res = falpha->GetParameter("Sigma");
0200 malpha_res = falpha->GetParameter("Mean");
0201
0202 cout << " " << wheel << " " << station << " " << sector << " salpha_res: " << salpha_res << " malpha_res: " << malpha_res;
0203 cout << " sbeta_res: " << sbeta_res << " mbeta_res: " << mbeta_res;
0204 cout << endl;
0205
0206 }
0207
0208 // Segment angle pull
0209 if (doPullabxy) {
0210 if (station!=4) {
0211 TH1F* tmpBeta = (TH1F*) hRes4D->hPullBetaRZ->Clone("tmpBeta");
0212 tmpBeta->Rebin(2);
0213 TF1* fbeta = drawGFit(tmpBeta, nsigma, -2. , 2.);
0214 sbeta_pull = fbeta->GetParameter("Sigma");
0215 mbeta_pull = fbeta->GetParameter("Mean");
0216
0217 TH1F* tmpY = (TH1F*) hRes4D->hPullYRZ->Clone("tmpY");
0218 tmpY->Rebin(2);
0219 TF1* fy = drawGFit(tmpY, nsigma, -2. , 2. );
0220 sy_pull = fy->GetParameter("Sigma");
0221 my_pull = fy->GetParameter("Mean");
0222 }
0223
0224 TH1F* tmpAlpha = (TH1F*) hRes4D->hPullAlpha->Clone("tmpAlpha");
0225 tmpAlpha->Rebin(2);
0226 TF1* falpha = drawGFit(tmpAlpha, nsigma, -2. , 2. );
0227 salpha_pull = falpha->GetParameter("Sigma");
0228 malpha_pull = falpha->GetParameter("Mean");
0229
0230 TH1F* tmpX = (TH1F*) hRes4D->hPullX->Clone("tmpX");
0231 tmpX->Rebin(2);
0232 TF1* fx = drawGFit(tmpX, nsigma, -2. , 2.);
0233 sx_pull = fx->GetParameter("Sigma");
0234 mx_pull = fx->GetParameter("Mean");
0235
0236
0237 cout << " " << wheel << " " << station << " " << sector << " salpha_pull: " << salpha_pull << " malpha_pull: " << malpha_pull;
0238 cout << " sbeta_pull: " << sbeta_pull << " mbeta_pull: " << mbeta_pull;
0239 cout << " sx_pull: " << sx_pull << " mx_pull: " << mx_pull;
0240 cout << " sy_pull: " << sy_pull << " my_pull: " << my_pull;
0241 cout << endl;
0242
0243 }
0244
0245
0246 double ratio=0;
0247
0248 if (doNSeg) {
0249 TH1F* hNs = hEff4D->hNSeg;
0250 double int_1 = hNs->Integral(2,21);
0251 double int_2 = hNs->Integral(3,21);
0252 ratio = int_2/int_1;
0253 // cout << "int_1: " << int_1 <<" int_2: " << int_2 << " ratio: " << ratio << endl;
0254 }
0255
0256
0257
0258 float p0_Phi=0;
0259 float p1_Phi=0;
0260 float p2_Phi=0;
0261 float p0_Theta=0;
0262 float p1_Theta=0;
0263 float p2_Theta=0;
0264
0265 if(doPar){
0266 float min_Phi;
0267 float max_Phi;
0268 float min_Theta;
0269 float max_Theta;
0270
0271 rangeAngle(abs(wheel), station, 1, min_Phi, max_Phi);
0272 TF1 *angleDep_Phi= new TF1("angleDep_Phi","[0]*cos(x+[2])+[1]", min_Phi, max_Phi);
0273 angleDep_Phi->SetParameters(0.01,0.001, 0.1);
0274
0275 TH1F* hprof_Phi = plotAndProfileX(hResPhi1->hResVsAngle,1,1,1,-0.04, 0.04, min_Phi-0.03, max_Phi+0.03);
0276 if((station!=1 ||(station==1 && wheel==0))) angleDep_Phi->FixParameter(2,0.);
0277 if((abs(wheel) == 1 ||wheel == -2) && station == 1) angleDep_Phi->SetParameter(2,-0.12);
0278
0279 hprof_Phi->Fit(angleDep_Phi,"RQN");
0280
0281 angleDep_Phi->Draw("same");
0282 p0_Phi = angleDep_Phi->GetParameter(0);
0283 p1_Phi = angleDep_Phi->GetParameter(1);
0284 p2_Phi = angleDep_Phi->GetParameter(2);
0285
0286 if (station!=4){
0287 rangeAngle(abs(wheel), station, 2, min_Theta, max_Theta);
0288 TF1 *angleDep_Theta= new TF1("angleDep_Theta","[0]*cos(x+[2])+[1]", min_Theta, max_Theta);
0289 angleDep_Theta->SetParameters(0.01,0.001, 0.1);
0290 angleDep_Theta->FixParameter(2,0.);
0291 TH1F* hprof_Theta = plotAndProfileX(hResTheta->hResVsAngle,1,1,1,-0.04, 0.04, min_Theta-0.03, max_Theta+0.03);
0292
0293 if((wheel==-2 || wheel==2) &&(station == 1 || station == 2)) {
0294 angleDep_Theta->FixParameter(0,0.);
0295 angleDep_Theta->FixParameter(1,0.);
0296 }
0297
0298 hprof_Theta->Fit(angleDep_Theta,"RQN");
0299 angleDep_Theta->Draw("same");
0300 //angleDep_Theta->Draw("same");
0301 p0_Theta = angleDep_Theta->GetParameter(0);
0302 p1_Theta = angleDep_Theta->GetParameter(1);
0303 p2_Theta = angleDep_Theta->GetParameter(2);
0304
0305 }
0306 }
0307
0308
0309 // Write summary table (one row per SL)
0310 int secmin=sector;
0311 int secmax=sector;
0312 if (sector ==0) {
0313 secmin=1;
0314 secmax=14;
0315 }
0316 for (int sec = secmin; sec<=secmax; sec++) {
0317 if (station!=4 && sec>12) continue;
0318 f << wheel << " " << station << " " << sec << " " << 1 << " " << effS1RPhi << " " << effS3RPhi << " " << effSeg << " " << sphi << " " << pphi << " " << malpha_res << " " <<salpha_res << " " << malpha_pull << " " <<salpha_pull << " " << mx_pull << " " <<sx_pull << " " << mx << " " << sx << " " <<ratio << " " << p0_Phi << " " << p1_Phi << " " << p2_Phi << endl;
0319 if (station!=4) f << wheel << " " << station << " " << sec << " " << 2 << " " << effS1RZ << " " << effS3RZ << " " << effSeg <<" " << stheta << " " << ptheta << " " << mbeta_res << " " <<sbeta_res << " " << mbeta_pull << " " <<sbeta_pull << " " << my_pull << " " <<sy_pull << " " << my << " " << sy << " " <<ratio <<" " << p0_Theta << " " << p1_Theta << " " << p2_Theta <<endl;
0320 f << wheel << " " << station << " " << sec << " " << 3 << " " << effS1RPhi << " " << effS3RPhi << " " << effSeg << " " << sphi << " " << pphi << " " << malpha_res << " " <<salpha_res << " " << malpha_pull << " " <<salpha_pull << " " << mx_pull << " " <<sx_pull << " " << mx << " " << sx << " " <<ratio << " " << p0_Phi << " " << p1_Phi << " " << p2_Phi <<endl;
0321 }
0322 } // sector
0323 } //station
0324 } //wheel
0325 }
