EcalSimAlgos/test/dumpProfile.C

0001 void dumpProfile() {
0002
0003  std::string fileName_ = "Profile_SM10.root";
0004  TFile *shapeFile_ = TFile::Open(fileName_.c_str(),"old");
0005  TProfile* PROF_704 = (TProfile*) shapeFile_->Get("SHAPE_XTAL_704");
0006
0007  ofstream out;
0008  out.open("dat.txt");
0009
0010  int nBinsHisto_ = 250;
0011  std::vector<double> shapeArray(nBinsHisto_,0.0);
0012
0013  double max = -999;
0014  int imax = 0;
0015  for(int ibin=0; ibin < nBinsHisto_; ++ibin)
0016    {
0017      out << "shapeArray[" << ibin << "] = " << PROF_704->GetBinContent(ibin+1) << " ; \n";
0018      shapeArray[ibin] = PROF_704->GetBinContent(ibin);
0019      std::cout << "Original shape, ns = " << ibin << " shape = " << shapeArray[ibin] << std::endl;
0020      if ( shapeArray[ibin] > max ) {
0021        max = shapeArray[ibin];
0022        imax = ibin;
0023      }
0024
0025    }//loop
0026
0027  out.close();
0028
0029  double xMinHisto_ = -1.;
0030  double xMaxHisto_ = 9.;
0031  double binw = (xMaxHisto_ - xMinHisto_)/(shapeArray.size());
0032  int nbins = shapeArray.size()/10;
0033
0034  float low =  xMinHisto_+(double)(imax-nbins/2+0.5)*binw;
0035  float up = xMinHisto_+(double)(imax+nbins/2+0.5)*binw;
0036
0037  double* x = new double[nbins];
0038  double* y = new double[nbins];
0039  for (int i = 0; i < nbins; i++) {
0040    x[i] = xMinHisto_ + (double)(imax - nbins/2 + i + 0.5)*binw;
0041    y[i] = shapeArray[imax - nbins/2 + i];
0042    std::cout << " x,y = " << x[i] << " " << y[i] << " " << (double)(imax - nbins/2 + i + 0.5) << std::endl;
0043  }
0044  TGraph* graph = new TGraph(nbins, x, y);
0045  graph->Fit("pol3", "V");//"Q 0");
0046  TF1* fFit = graph->GetFunction("pol3");
0047  double tMax = fFit->GetMaximumX();
0048
0049  std:;cout << "Maxiumum = " << tMax << std::endl;
0050
0051  gStyle->SetOptFit(1111);
0052
0053  TCanvas *MyC = new TCanvas("MyC","Test canvas",1);
0054  MyC->Divide(2,1);
0055  MyC->cd(1);
0056  PROF_704->Draw();
0057  MyC->cd(2);
0058  fFit->Draw();
0059  MyC->SaveAs("PROF_704.jpg");
0060
0061 }