#include "TPaveStats.h"
#include "TCanvas.h"
#include "TStyle.h"
#include "TLatex.h"
#include "TH2F.h"
#include "TCanvas.h"
#include "TFile.h"
#include "TObject.h"
#include "TGaxis.h"
#include "TAxis.h"
#include "TPaletteAxis.h"
#include "TPad.h"
#include "TPolyLine.h"
#include "TASImage.h"
#include <iostream>
#include <fstream>
using namespace std;
void AnalysePSFReflector()
{
//Prepare some global variables 
  
  char line[16], file[4], pos[4];
  int logic,posx,posy,x,y,maxx,maxy,maxz,maxzglobal=-1000,minz,minzglobal=1000,counter=0;
  float distance[17][17],r90[17][17],r80[17][17],r50[17][17],integralr50[17][17],integralr80[17][17],integralr90[17][17],r90max=0.1,r90min=40,r90mean=0,r80max=0.1,r80min=40,r80mean=0.,r50max,r50min,r50mean=0.,integralr50max,integralr50min,integralr50mean,integralr80max,integralr80min,integralr80mean,integralr90max,integralr90min,integralr90mean; 
 //Prepare positions on canvas referring to mirror position at reflector
  int position[17][17]={145,128,111,94,77,60,43,26,9,162,179,196,213,230,247,264,281,
			146,129,112,95,78,61,44,27,10,163,180,197,214,231,248,265,282,
			147,130,113,96,79,62,45,28,11,164,181,198,215,232,249,266,283,
			148,131,114,97,80,63,46,29,12,165,182,199,216,233,250,267,284,
			149,132,115,98,81,64,47,30,13,166,183,200,217,234,251,268,285,
			150,133,116,99,82,65,48,31,14,167,184,201,218,235,252,269,286,
			151,134,117,100,83,66,49,32,15,168,185,202,219,236,253,270,287,
			152,135,118,101,84,67,50,33,16,169,186,203,220,237,254,271,288,
			153,136,119,102,85,68,51,34,17,170,187,204,221,238,255,272,289,
			144,127,110,93,76,59,42,25,8,161,178,195,212,229,246,263,280,
			143,126,109,92,75,58,41,24,7,160,177,194,211,228,245,262,279,
			142,125,108,91,74,57,40,23,6,159,176,193,210,227,244,261,278,
			141,124,107,90,73,56,39,22,5,158,175,192,209,226,243,260,277,
			140,123,106,89,72,55,38,21,4,157,174,191,208,225,242,259,276,
			139,122,105,88,71,54,37,20,3,156,173,190,207,224,241,258,275,
			138,121,104,87,70,53,36,19,2,155,172,189,206,223,240,257,274,
			137,120,103,86,69,52,35,18,1,154,171,188,205,222,239,256,273};
  //Prepare polyline along pad's XlowNDC and YlowNDC coordinates in order to distinguish between aluminium and galss mirrors//pad->cd(padnumber)->GetXlowNDC() and GetYlowNDC() delivered the coordinates written in the Double_t arrays below
  Double_t plinex[29] = { 0.296689,0.703817,0.703817,0.761979,0.761979,0.82014,0.82014,0.878301,0.878301,0.82014,0.82014,0.761979,0.761979,0.703817,0.703817,0.296689,0.296689,0.238527,0.238527,0.180366,0.180366,0.122205,0.122205,0.180366,0.180366,0.238527,0.238527,0.296689,0.296689};
  Double_t pliney[29] = { 0.877852,0.877852,0.819721,0.819721,0.76159,0.76159,0.703459,0.703459,0.296541,0.296541,0.23841,0.23841,0.180279,0.180279,0.122148,0.122148,0.180279,0.180279,0.23841,0.23841,0.296541,0.296541,0.703459,0.703459,0.76159,0.76159,0.819721,0.819721,0.877852};
  TPolyLine *pline = new TPolyLine(29,plinex,pliney);
  pline->SetLineWidth(2);
  //Prepare file list into which al existing data bla are written
  listFile=string("Files.txt");
  ifstream ffile(listFile.c_str() );
  if( !ffile ) {
    cout <<"Error opening " << listFile << " for input. Creating it!" << endl;
    //create file list 
    system("ls mirror????.root > Files.txt");
  }
    //Prepare canvas
    TCanvas *c0 = new TCanvas("c0","PSF of MAGIC Reflector",1000,1000);
    gStyle->SetOptStat(0);
    gStyle->SetOptTitle(0);
    c0->SetFillColor(0);
    c0->SetBorderSize(0);
    //Add pads where magic reflector, z scale and legend are shown, p4 is divided so that data can be arranged as mirrors are installed in reflector
    TPad *p2= new TPad("p2","Statistics",0,0.77,0.11,0.93);
    TPad *p3= new TPad("p3","Z Scale",0.95,0,1,0.93);
    TPad *p4= new TPad("p4","Reflector Scheme",0,0,0.95,0.93);
    TPad *paddummy= new TPad("paddummy","Reflector Scheme",0,0,0.95,0.93);
    p3->Draw();
    p4->Divide(17,17,0,0);
    p4->Draw();
    paddummy->Draw();
    paddummy->SetFillStyle(4000);  
    // Draw picture title
    TPaveLabel *title1 = new TPaveLabel(0.06,0.95,0.86,0.99,"PSF of the MAGIC II mirrors");
    title1->SetFillColor(0);
    title1->SetBorderSize(0);
    title1->SetTextFont(52);
    title1->Draw();
    //Open file list
    FILE *fpdata=fopen("Files.txt","r");
    //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    //Draw PSF in graphs
      while (!feof(fpdata)){
	fscanf(fpdata,"%s",&line);
	//Get mirror position from .root file name
	strncpy(pos,&line[6],4);
	//Set Ende-sign manually
	pos[4] = '\0';
	logic=strcmp(pos,file);
	if(logic!=0){
	  strncpy(file,&pos[0],4);
	  file[4] = '\0';
	  //Separte positions to x and y coordinate
	  sscanf(pos,"%d %d",&posx,&posy);
	  x=posx;
	  y=posy;
	  //Convert negative coordinates in order to retrieve position from the array
	  if(x<0)
	    {
	      x=TMath::Abs(x)+8;
	    }
	  if(y<0)
	    {
	      y=TMath::Abs(y)+8;
	    } 
	  //file whether data file does really exist
	  FILE *fp;
	  if ((fp=fopen(line,"r"))!=NULL) {
	    //Open data file as .root file
	    TString *string1=new TString(line,15);
	    TFile *tfile= TFile::Open(*string1,"R");
	    //Retrieve single suppads from canvas
	    TCanvas *c1 = ( TCanvas*)tfile->Get("Data presentation");
	    int width=c1->GetWw();
	    TList *prims=(TList*)c1->GetListOfPrimitives();
	    int nprims = prims->GetSize();
	    //Some advice by root support how to be able to delete pointers and data structure after use successfully by becoming owner of primitives of each pad
	    for (int i=0;i<nprims;i++) {
	      if (prims->At(i)->InheritsFrom(TPad::Class())) {
		TPad *pad = (TPad*)prims->At(i);
		pad->GetListOfPrimitives()->SetOwner();
	      }
	    }
	    //Get graph as well as R50,R80,R90 as well as integrals from pad by using x position of indicating line
	    TPad *p6 =(TPad*)c1->FindObject("Data presentation_9");
	    TList *list2=(TList*)p6->GetListOfPrimitives();
	    TGraph *g0=(TGraph*)list2->FindObject("");
	    TLine *o2=(TLine*)list2->FindObject("TLine");
	    TLine *o3=(TLine*)list2->At(5);
	    TLine *o4=(TLine*)list2->At(8);
	    //Select whether r50, r90 and r80 have a reasonable value Calculate max and min of r50, r90 and r80 as well as average, calculate respective integrals//Due to data analysis procedure change kind of confused computing //STUPID CORNELIA for later results sufficient to use GetX1 method	   
	    r50[x][y]=o2->GetX1();
	    r90[x][y]=o3->GetX1();
	    r80[x][y]=o4->GetX1();
	    integralr50[x][y]=o2->GetY2();
	    integralr90[x][y]=o3->GetY2();
	    integralr80[x][y]=o4->GetY2();
	    if(r90max<r90[x][y]){
	      r90max=r90[x][y];
	    }
	    if(r90min>r90[x][y]){
	      r90min=r90[x][y];
	    }
	    if(r80max<r80[x][y]){
	      r80max=r80[x][y];
	    }
	    if(r80min>r80[x][y]){
	      r80min=r80[x][y];
	    }
	    if(r50max<r50[x][y]){
	      r50max=r50[x][y];
	    }
	    if(r50min>r50[x][y]){
	      r50min=r50[x][y];
	    }
	    if(integralr90max<integralr90[x][y]){
	      integralr90max=integralr90[x][y];
	    }
	    if(integralr90min>integralr90[x][y]){
	      integralr90min=integralr90[x][y];  
	    }
	    if(integralr80max<integralr80[x][y]){
	      integralr80max=integralr80[x][y];
	    }
	    if(integralr80min>integralr80[x][y]){
	      integralr80min=integralr80[x][y];  
	    }
	    if(integralr50max<integralr50[x][y]){
	      integralr50max=integralr50[x][y];
	    }
	    if(integralr50min>integralr50[x][y]){
	      integralr50min=integralr50[x][y];
	    }
	    counter++;
	    r90mean+=r90[x][y];
	    r80mean+=r80[x][y];
	    r50mean+=r80[x][y];
	    integralr90mean+=integralr90[x][y];
	    integralr80mean+=integralr80[x][y];
	    integralr50mean+=integralr50[x][y];
	    p4->cd(position[x][y]);
	    //p4->cd(position[x][y])->SetFrameBorderMode(1);
	    g0->GetXaxis()->SetTitle("");
	    g0->GetXaxis()->SetLabelSize(0);
	    g0->GetXaxis()->SetNdivisions(0);
	    g0->GetYaxis()->SetTitleSize(0);
	    g0->GetYaxis()->SetLabelSize(0);
	    g0->GetYaxis()->SetNdivisions(0);
	    g0->DrawClone("AP");
	    p4->cd(position[x][y])->Update();
	    delete c1;
	    delete string1;
	    delete tfile;
	    fclose(fp);
	  }
	}
      }
      r90mean=r90mean/counter;
      r80mean=r80mean/counter;
      r50mean=r50mean/counter;
      integralr90mean=integralr90mean/counter;
      integralr80mean=integralr80mean/counter;
      integralr50mean=integralr50mean/counter;
      p3->cd();
      char text[600];
      sprintf(text,"#splitline{#splitline{#splitline{#splitline{#splitline{R_{90Max}:%.1f mm}{R_{80Max}:%.1f mm}}{#splitline{R_{50Max}:%.1f mm}{R_{90Min}:%.1f mm}}}{#splitline{#splitline{R_{80Min}:%.1f mm}{R_{50Min}:%.1f mm}}{#splitline{R_{90Mean}:%.1f mm}{R_{90Mean}:%.1f mm}}}}{#splitline{#splitline{#splitline{R_{50Men}:%.1f mm}{I_{R_{90Max}}:%.3e}}{#splitline{I_{R_{90Max}}:%.3e}{I_{R_{50Min}}:%.3e}}}{#splitline{#splitline{I_{R_{90Min}}:%.3e}{I_{R_{80Min}}:%.3e}}{#splitline{I_{R_{50Min}}:%.3e}{I_{R_{90Mean}}:%.3e}}}}}{#splitline{I_{R_{80Mean}}:%.3e}{I_{R_{50Mean}}:%.3e}}",r90max,r80max,r50max,r90min,r80min,r50min,r90mean, r80mean, r50mean,integralr90max,integralr80max,integralr50max,integralr90min,integralr80min,integralr50min,integralr90mean,integralr80mean,integralr50mean);
      TLatex *text1=new TLatex(0.01,0.33, text);
      text1->SetTextSize(0.14);
      text1->Draw();
      //Draw polyline
      paddummy->cd();
      pline->SetLineColor(2);
      pline->Draw(); 
       char save[50];
   
      sprintf(save,"PSFinGraphs.ps");
      c0->Print(save);
      delete title1;
    delete p2;
    delete p3;
    delete p4;
    delete c0;
  }
}