EMCal/EMCALmacros/ProjectionReverse2.C

0001 #include <cmath>
0002 #include <TFile.h>
0003 #include <TString.h>
0004 #include <TLine.h>
0005 #include <TTree.h>
0006 #include <TLatex.h>
0007 #include <TGraphErrors.h>
0008 #include <cassert>
0009 #include <iostream>
0010 #include <fstream>
0011 using namespace std;
0012
0013
0014
0015 using std::cout;
0016 using std::endl;
0017 #endif
0018
0019 void ProjectionReverse2()
0020 {
0021
0022      gStyle->SetOptFit(0);
0023    gStyle->SetOptStat(0);
0024
0025     TFile *fin = new TFile("EnergyPosition2.root");
0026
0027 double y1 = -90;
0028 double y2 = -117;
0029 double y3 = -143;
0030 double y4 = -167;
0031 double y5 = -192;
0032
0033 //Set Tower Vertical Lines//
0034
0035 double x1 = -205;
0036 double x2 = -228;
0037 double x3 = -253;
0038 double x4 = -277;
0039 double x5 = -301;
0040
0041
0042
0043 int binx;
0044
0045 int forward = 2;
0046 int backward = 2;
0047
0048
0049 int biny;
0050 double step = 1;
0051 /*
0052 double Xmin = 170;
0053 double Xmax = 340.0;
0054
0055 double Ymin = 60;
0056 double Ymax = 245.0;
0057 */
0058 double Xmax = -170;
0059 double Xmin = -310.0;
0060
0061 double Ymax = -60;
0062 double Ymin = -205.0;
0063 cout << Ymax << endl;
0064
0065 int NumX = (Xmax - Xmin)/step;
0066 int NumY = (Ymax - Ymin)/step;
0067
0068 const int XBins = NumX;
0069
0070 const int YBins = NumY;
0071 double x;
0072 double y;
0073
0074 cout << "XBins = " << XBins << endl;
0075 cout << "YBins = " << YBins << endl;
0076
0077
0078 double value[XBins][YBins];
0079
0080 double finalvalue[XBins][YBins];
0081
0082 TH2D *Inter = new TH2D("Inter","",XBins,Xmin,Xmax,YBins,Ymin,Ymax);
0083
0084 /*
0085 double ylower = 35;
0086 double yhigher = 188;
0087 double xlower = 178;
0088 double xhigher = 301;
0089 */
0090 int k;
0091
0092 int NumberofValue;
0093 double Sum;
0094 int Xnow;
0095
0096 int Ynow;
0097 TCanvas *c1 = new TCanvas("c1", "c1",0,0,800,600);
0098
0099 c1->cd();
0100
0101 EnPo->GetXaxis()->SetTitle("Horizontal Axis (mm)");
0102 EnPo->GetYaxis()->SetTitle("Vertical Axis (mm)");
0103
0104 EnPo->SetTitle("Energy vs Horizontal and Vertical Positions Before Interpolation - 0 Degree");
0105
0106 EnPo->Draw("colz");
0107 c1->Update();
0108
0109 c1->SaveAs("Before0Degree.png");
0110
0111 for(int i = 0; i < XBins; i++)
0112 {
0113
0114 //x = Xmin + step*i;
0115
0116
0117 for(int j = 0; j < YBins; j++)
0118 {
0119
0120 //  TCanvas *c1 = new TCanvas("c1", "c1",0,0,800,600);
0121 //  y = Ymin + step*j;
0122
0123
0124
0125
0126 value[i][j] = EnPo->GetBinContent(i,j);
0127
0128
0129 }
0130 }
0131
0132 //Horizontal Interpolation Begins//
0133 for(int i = backward; i < XBins - backward; i++)
0134 {
0135
0136     for(int j = backward; j < YBins - backward; j++)
0137 {
0138
0139 //cout << "OK" << endl;
0140
0141 NumberofValue = 0;
0142 Sum = 0;
0143
0144 Xnow = XBins - i;
0145 Ynow = YBins - j;
0146
0147
0148 for(int k = i - backward; k < i + forward +1; k++)
0149 {
0150
0151     for(int m = j - backward; m < j+backward+1; m++)
0152     {
0153
0154         if(value[k][m] > 0) NumberofValue = NumberofValue + 1;
0155
0156         Sum = Sum + value[k][m];
0157
0158             }
0159
0160 }
0161
0162 cout << "Number of Value = " << NumberofValue << endl;
0163
0164 cout << "Sum = " << Sum << endl;
0165
0166 if(NumberofValue > 0) finalvalue[i][j] = Sum/NumberofValue;
0167
0168 cout << " i = " << i <<  "  j = " << j << endl;
0169 cout << " Xnow = " << Xnow <<  "  YNow = " << Ynow << endl;
0170
0171 cout << "Final Value is " << finalvalue[i][j]  << endl;
0172
0173
0174 //Inter->SetBinContent(i,j,finalvalue[i][j]);
0175
0176 Inter->SetBinContent(Xnow,Ynow,finalvalue[i][j]);
0177
0178 }
0179
0180
0181 }
0182
0183 /*
0184 cout << "Horizontal Done" << endl;
0185
0186
0187 //Vertical Interpolation Begins//
0188
0189 cout << "DONE Value" << endl;
0190
0191 for(int i = 1; i < XBins; i++)
0192 {
0193
0194     for(int j = 1; j < YBins; j++)
0195 {
0196
0197     cout << "Before Updating: Binx =  " << i << "   " << "Bin y = " << j << " Value =  " << value[i][j] << endl;
0198
0199
0200 if (value[i][j] < 6.8 && value[i][j-1] >0 && i < XBins && j < YBins)
0201 {
0202
0203
0204 k = j;
0205 do {
0206     k = k + 1;
0207
0208 //  cout << "k value is = " << k << endl;
0209 }
0210 while(value[i][k] == 0 && k < YBins);
0211
0212 //cout << "neW k " << k  << endl;
0213
0214 value[i][j] =value[i][j-1] + (value[i][k] - value[i][j-1])/(k+1-j);
0215
0216     cout << "After Updating: Binx =  " << i << "   " << "  Bin y = " << j << "  Value =  " <<value[i][j] << endl;
0217
0218 }
0219
0220
0221 }
0222
0223
0224 }
0225
0226 cout << "Vk =" << k << endl;
0227
0228 //Vertical Interpolation Ends//
0229
0230 cout << "VERTICAL DONE" << endl;
0231
0232 */
0233 /*
0234 for(int i =0; i < XBins; i++)
0235 {
0236
0237
0238
0239     for(int j = 1; j < YBins; j++)
0240 {
0241
0242
0243 }
0244
0245
0246 }
0247 */
0248 Inter->GetXaxis()->SetTitle("Horizontal Axis (mm)");
0249 Inter->GetYaxis()->SetTitle("Vertical Axis (mm)");
0250 Inter->SetTitle("Energy vs Horizontal and Vertical Positions After Interpolation - 10 Degree");
0251
0252   TFile *fout = new TFile("Interpolated2.root","RECREATE");
0253 Inter->Write();
0254
0255
0256
0257
0258
0259
0260
0261 Inter->Draw("colz");
0262
0263
0264 char Outname[512];
0265
0266 const int All = 7;
0267 double inteval = -12.5;
0268 /*
0269 double left;
0270 double right;
0271 double initial = -191;
0272
0273
0274 for(int i=0; i < All; i++)
0275 {
0276 TCanvas *c22 = new TCanvas("c22", "c22",0,0,800,600);
0277
0278 Inter->SetMinimum(3.5);
0279 Inter->SetMaximum(9.0);
0280 c22->cd();
0281     Inter->Draw("colz");
0282
0283     left  = initial + inteval * i - 5;
0284     right  = initial + inteval * i + 5;
0285
0286    TLine *l1 = new TLine(left,Ymin,left,Ymax);
0287
0288    l1->SetLineWidth(2);
0289    l1->Draw("same");
0290
0291    TLine *l2 = new TLine(right,Ymin,right,Ymax);
0292    l2->SetLineWidth(2);
0293    l2->Draw("same");
0294 sprintf(Outname,"Result/VerticalPlot%d.png",i);
0295
0296 c22->Update();
0297
0298 c22->SaveAs(Outname);
0299
0300 }
0301 */
0302 double initial = -77;
0303 int HAll = All + 5;
0304 for(int i=0; i < HAll; i++)
0305 {
0306 TCanvas *c22 = new TCanvas("c22", "c22",0,0,800,600);
0307
0308 Inter->SetMinimum(3.5);
0309 Inter->SetMaximum(9.0);
0310 c22->cd();
0311     Inter->Draw("colz");
0312
0313     left  = initial + inteval * i - 5;
0314     right  = initial + inteval * i + 5;
0315
0316    TLine *l1 = new TLine(Xmin,left,Xmax,left);
0317
0318    l1->SetLineWidth(2);
0319    l1->Draw("same");
0320
0321    TLine *l2 = new TLine(Xmin,right,Xmax,right);
0322    l2->SetLineWidth(2);
0323    l2->Draw("same");
0324 sprintf(Outname,"Result2/HorizontalPlot%d.png",i);
0325
0326 c22->Update();
0327
0328 c22->SaveAs(Outname);
0329
0330 }
0331
0332
0333 /*
0334 // Horizontal line Tower 1 //
0335    TLine *l1 = new TLine(Xmin,y1,Xmax,y1);
0336    l1->SetLineWidth(1);
0337    l1->Draw("same");
0338
0339
0340  TLatex* texCol;
0341 texCol= new TLatex(Xmin,y1,"Tower 1");
0342   texCol->SetTextAlign(13);
0343   texCol->SetTextSize(0.03);
0344  texCol->SetTextFont(42);
0345   texCol->Draw();
0346
0347
0348 // Horizontal line Tower 2 //
0349
0350
0351    TLine *l2 = new TLine(Xmin,y2,Xmax,y2);
0352    l2->SetLineWidth(1);
0353    l2->Draw("same");
0354
0355
0356
0357 texCol= new TLatex(Xmin,y2,"Tower 2");
0358   texCol->SetTextAlign(13);
0359   texCol->SetTextSize(0.03);
0360  texCol->SetTextFont(42);
0361   texCol->Draw();
0362
0363 // Horizontal line Tower 3//
0364
0365
0366    TLine *l3 = new TLine(Xmin,y3,Xmax,y3);
0367    l3->SetLineWidth(1);
0368    l3->Draw("same");
0369
0370
0371
0372 texCol= new TLatex(Xmin,y3,"Tower 3");
0373   texCol->SetTextAlign(13);
0374   texCol->SetTextSize(0.03);
0375  texCol->SetTextFont(42);
0376   texCol->Draw();
0377
0378
0379
0380 // Horizontal line Tower 4 //
0381
0382
0383    TLine *l4 = new TLine(Xmin,y4,Xmax,y4);
0384    l4->SetLineWidth(1);
0385    l4->Draw("same");
0386
0387
0388
0389 texCol= new TLatex(Xmin,y4,"Tower 4");
0390   texCol->SetTextAlign(13);
0391   texCol->SetTextSize(0.03);
0392  texCol->SetTextFont(42);
0393   texCol->Draw();
0394
0395   // Horizontal line Tower 5 //
0396
0397
0398    TLine *l4 = new TLine(Xmin,y5,Xmax,y5);
0399    l4->SetLineWidth(1);
0400    l4->Draw("same");
0401
0402
0403
0404
0405 texCol= new TLatex(Xmin,y5,"Tower 5");
0406   texCol->SetTextAlign(13);
0407   texCol->SetTextSize(0.03);
0408  texCol->SetTextFont(42);
0409   texCol->Draw();
0410
0411
0412
0413
0414 // Vertical line Tower 1 //
0415    TLine *l6 = new TLine(x1,Ymin,x1,Ymax);
0416    l6->SetLineWidth(1);
0417    l6->Draw("same");
0418
0419
0420  TLatex* texCol;
0421 texCol= new TLatex(x1,Ymax,"Tower 1");
0422   texCol->SetTextAlign(33);
0423   texCol->SetTextSize(0.03);
0424  texCol->SetTextFont(42);
0425   texCol->Draw();
0426
0427
0428 // Horizontal line Tower 2 //
0429
0430
0431    TLine *l7 = new TLine(x2,Ymin,x2,Ymax);
0432    l7->SetLineWidth(1);
0433    l7->Draw("same");
0434
0435
0436  TLatex* texCol;
0437 texCol= new TLatex(x2,Ymax,"Tower 2");
0438   texCol->SetTextAlign(33);
0439   texCol->SetTextSize(0.03);
0440  texCol->SetTextFont(42);
0441   texCol->Draw();
0442
0443 // Horizontal line Tower 3//
0444
0445
0446    TLine *l8 = new TLine(x3,Ymin,x3,Ymax);
0447    l8->SetLineWidth(1);
0448    l8->Draw("same");
0449
0450
0451  TLatex* texCol;
0452 texCol= new TLatex(x3,Ymax,"Tower 3");
0453   texCol->SetTextAlign(33);
0454   texCol->SetTextSize(0.03);
0455  texCol->SetTextFont(42);
0456   texCol->Draw();
0457
0458
0459 // Horizontal line Tower 4 //
0460
0461    TLine *l9 = new TLine(x4,Ymin,x4,Ymax);
0462    l9->SetLineWidth(1);
0463    l9->Draw("same");
0464
0465
0466  TLatex* texCol;
0467 texCol= new TLatex(x4,Ymax,"Tower 4");
0468   texCol->SetTextAlign(33);
0469   texCol->SetTextSize(0.03);
0470  texCol->SetTextFont(42);
0471   texCol->Draw();
0472
0473
0474   // Horizontal line Tower 5 //
0475
0476
0477    TLine *l10 = new TLine(x5,Ymin,x5,Ymax);
0478    l10->SetLineWidth(1);
0479    l10->Draw("same");
0480
0481
0482  TLatex* texCol;
0483 texCol= new TLatex(x5,Ymax,"Tower 5");
0484   texCol->SetTextAlign(33);
0485   texCol->SetTextSize(0.03);
0486  texCol->SetTextFont(42);
0487   texCol->Draw();
0488 */
0489
0490
0491
0492
0493
0494
0495
0496
0497 }