EMCal/EMCALmacros/ProjectionReverse.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 ProjectionReverse()
0020 {
0021
0022      gStyle->SetOptFit(0);
0023    gStyle->SetOptStat(0);
0024
0025     TFile *fin = new TFile("EnergyPosition3.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 = -340.0;
0060
0061 double Ymax = -60;
0062 double Ymin = -245.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 - 0 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
0268 double left;
0269 double right;
0270 double initial = -191;
0271 double inteval = -12.5;
0272
0273 for(int i=0; i < All; i++)
0274 {
0275 TCanvas *c22 = new TCanvas("c22", "c22",0,0,800,600);
0276
0277 Inter->SetMinimum(3.5);
0278 Inter->SetMaximum(9.0);
0279 c22->cd();
0280     Inter->Draw("colz");
0281
0282     left  = initial + inteval * i - 5;
0283     right  = initial + inteval * i + 5;
0284
0285    TLine *l1 = new TLine(left,Ymin,left,Ymax);
0286
0287    l1->SetLineWidth(2);
0288    l1->Draw("same");
0289
0290    TLine *l2 = new TLine(right,Ymin,right,Ymax);
0291    l2->SetLineWidth(2);
0292    l2->Draw("same");
0293 sprintf(Outname,"Result/VerticalPlot%d.png",i);
0294
0295 c22->Update();
0296
0297 c22->SaveAs(Outname);
0298
0299 }
0300
0301 initial = -77;
0302 int HAll = All + 5;
0303 for(int i=0; i < HAll; i++)
0304 {
0305 TCanvas *c22 = new TCanvas("c22", "c22",0,0,800,600);
0306
0307 Inter->SetMinimum(3.5);
0308 Inter->SetMaximum(9.0);
0309 c22->cd();
0310     Inter->Draw("colz");
0311
0312     left  = initial + inteval * i - 5;
0313     right  = initial + inteval * i + 5;
0314
0315    TLine *l1 = new TLine(Xmin,left,Xmax,left);
0316
0317    l1->SetLineWidth(2);
0318    l1->Draw("same");
0319
0320    TLine *l2 = new TLine(Xmin,right,Xmax,right);
0321    l2->SetLineWidth(2);
0322    l2->Draw("same");
0323 sprintf(Outname,"Result/HorizontalPlot%d.png",i);
0324
0325 c22->Update();
0326
0327 c22->SaveAs(Outname);
0328
0329 }
0330
0331
0332 /*
0333 // Horizontal line Tower 1 //
0334    TLine *l1 = new TLine(Xmin,y1,Xmax,y1);
0335    l1->SetLineWidth(1);
0336    l1->Draw("same");
0337
0338
0339  TLatex* texCol;
0340 texCol= new TLatex(Xmin,y1,"Tower 1");
0341   texCol->SetTextAlign(13);
0342   texCol->SetTextSize(0.03);
0343  texCol->SetTextFont(42);
0344   texCol->Draw();
0345
0346
0347 // Horizontal line Tower 2 //
0348
0349
0350    TLine *l2 = new TLine(Xmin,y2,Xmax,y2);
0351    l2->SetLineWidth(1);
0352    l2->Draw("same");
0353
0354
0355
0356 texCol= new TLatex(Xmin,y2,"Tower 2");
0357   texCol->SetTextAlign(13);
0358   texCol->SetTextSize(0.03);
0359  texCol->SetTextFont(42);
0360   texCol->Draw();
0361
0362 // Horizontal line Tower 3//
0363
0364
0365    TLine *l3 = new TLine(Xmin,y3,Xmax,y3);
0366    l3->SetLineWidth(1);
0367    l3->Draw("same");
0368
0369
0370
0371 texCol= new TLatex(Xmin,y3,"Tower 3");
0372   texCol->SetTextAlign(13);
0373   texCol->SetTextSize(0.03);
0374  texCol->SetTextFont(42);
0375   texCol->Draw();
0376
0377
0378
0379 // Horizontal line Tower 4 //
0380
0381
0382    TLine *l4 = new TLine(Xmin,y4,Xmax,y4);
0383    l4->SetLineWidth(1);
0384    l4->Draw("same");
0385
0386
0387
0388 texCol= new TLatex(Xmin,y4,"Tower 4");
0389   texCol->SetTextAlign(13);
0390   texCol->SetTextSize(0.03);
0391  texCol->SetTextFont(42);
0392   texCol->Draw();
0393
0394   // Horizontal line Tower 5 //
0395
0396
0397    TLine *l4 = new TLine(Xmin,y5,Xmax,y5);
0398    l4->SetLineWidth(1);
0399    l4->Draw("same");
0400
0401
0402
0403
0404 texCol= new TLatex(Xmin,y5,"Tower 5");
0405   texCol->SetTextAlign(13);
0406   texCol->SetTextSize(0.03);
0407  texCol->SetTextFont(42);
0408   texCol->Draw();
0409
0410
0411
0412
0413 // Vertical line Tower 1 //
0414    TLine *l6 = new TLine(x1,Ymin,x1,Ymax);
0415    l6->SetLineWidth(1);
0416    l6->Draw("same");
0417
0418
0419  TLatex* texCol;
0420 texCol= new TLatex(x1,Ymax,"Tower 1");
0421   texCol->SetTextAlign(33);
0422   texCol->SetTextSize(0.03);
0423  texCol->SetTextFont(42);
0424   texCol->Draw();
0425
0426
0427 // Horizontal line Tower 2 //
0428
0429
0430    TLine *l7 = new TLine(x2,Ymin,x2,Ymax);
0431    l7->SetLineWidth(1);
0432    l7->Draw("same");
0433
0434
0435  TLatex* texCol;
0436 texCol= new TLatex(x2,Ymax,"Tower 2");
0437   texCol->SetTextAlign(33);
0438   texCol->SetTextSize(0.03);
0439  texCol->SetTextFont(42);
0440   texCol->Draw();
0441
0442 // Horizontal line Tower 3//
0443
0444
0445    TLine *l8 = new TLine(x3,Ymin,x3,Ymax);
0446    l8->SetLineWidth(1);
0447    l8->Draw("same");
0448
0449
0450  TLatex* texCol;
0451 texCol= new TLatex(x3,Ymax,"Tower 3");
0452   texCol->SetTextAlign(33);
0453   texCol->SetTextSize(0.03);
0454  texCol->SetTextFont(42);
0455   texCol->Draw();
0456
0457
0458 // Horizontal line Tower 4 //
0459
0460    TLine *l9 = new TLine(x4,Ymin,x4,Ymax);
0461    l9->SetLineWidth(1);
0462    l9->Draw("same");
0463
0464
0465  TLatex* texCol;
0466 texCol= new TLatex(x4,Ymax,"Tower 4");
0467   texCol->SetTextAlign(33);
0468   texCol->SetTextSize(0.03);
0469  texCol->SetTextFont(42);
0470   texCol->Draw();
0471
0472
0473   // Horizontal line Tower 5 //
0474
0475
0476    TLine *l10 = new TLine(x5,Ymin,x5,Ymax);
0477    l10->SetLineWidth(1);
0478    l10->Draw("same");
0479
0480
0481  TLatex* texCol;
0482 texCol= new TLatex(x5,Ymax,"Tower 5");
0483   texCol->SetTextAlign(33);
0484   texCol->SetTextSize(0.03);
0485  texCol->SetTextFont(42);
0486   texCol->Draw();
0487 */
0488
0489
0490
0491
0492
0493
0494
0495
0496 }