#include <stdio.h>

// if mode == 0 the macro exit after having printed the maps

void run()
{

  cout << "How many sets of file do you want to \"sum\"?" << endl;

  Int_t i;

  cin >> i;

  const Int_t dim = (const Int_t) i;

  cout << "Type the output path" << endl;

  std::string OutPath;

  cin >> OutPath;

  cout << "Type a root for the output names" << endl;

  std::string OutName;

  cin >> OutName;

  cout << "Type the title for the final plot" << endl;
  
  std::string FinalTitle;

  cin >> FinalTitle;

  std::string nomi[dim];
  std::string title[dim];

  // executing analysis on the selected files

  for (Int_t j=0 ; j < dim ; j++)
    {
      
      cout << j+1 << "^ file to be opened: " << endl;
      cin >> nomi[j];
      cout << "Title for this data set (no blank spaces): " << endl;
      cin >> title[j];
      cout << "---------------------------------------------------" << endl;

    }

    gROOT->Macro("$VMCWORKDIR/gconfig/rootlogon.C");
    
    TCanvas *can1 = new TCanvas();
    TCanvas *can2 = new TCanvas();

    Int_t max_bin,n_bin;

    Double_t step_size;

    max_bin = 1.;

    n_bin = 100000;

    step_size = max_bin / n_bin;

    TH3F *mapMat = new TH3F("Map","Map",500,-3.14,+3.14,500,-100,+60,2000,0,1.);
    //TH3F *mapMat = new TH3F("Map","Map",500,-3.14,+3.14,500,-100,+60,20000,0,3.);

    TH2F *mapMat2 = new TH2F("z scan","z scan",500,-100,+60,n_bin, -0.5 * step_size , max_bin - 0.5 * step_size);

    //    TH2F *mapMat2 = new TH2F("z scan","z scan",500,-100,+60,100000,0.,1.);

    //    TH2F *mapMat2 = new TH2F("z scan","z scan",500,-100,+60,2000,0.,6.);

    THStack *hs = new THStack("Stacked histos","test");

    leg = new TLegend(0.15,0.65,0.40,0.89);

    for (Int_t w = 0 ; w < dim ; w++)
      {

	mapMat2->Reset();

	    CombPlots(nomi[w],mapMat,mapMat2);

	    TProfile *pro = mapMat2->ProfileX();

	    TH1D *sum = pro->ProjectionX();

	    sum->SetFillColor(2+w);
	    leg->AddEntry(sum,(title[w]).c_str(),"F");
	    
	    hs->Add(sum);

      }

    
    // Plotting results

    // gStyle->SetPalette(1.);
    
    //    TProfile2D *proj = mapMat->Project3DProfile("yx");
    
    //    can1->cd();
    
    //    proj->Draw("COLZ");
    
    can2->cd();
    
      
    //sum->Draw("HIST");
    
    std::string figMap = Form("%sMap.png",OutName);
    std::string figProf = Form("%sProf.png",OutName);
    std::string figMapLOG = Form("%sMapLOG.png",OutName);
    
    //  can1->Print(figMap.c_str(),"png");
    //  can2->Print(figProf.c_str(),"png");
  
    can1->cd();
    
    //    gPad->SetLogz(1);
    
    //     can1->Print(figMapLOG.c_str(),"png");
    
    //    if (mode==0) exit();    


    hs->Draw("HIST");

    hs->GetXaxis()->SetTitle("Z [cm]");
    hs->GetYaxis()->SetTitle("# Rad Lengths");
    hs->SetTitle(FinalTitle.c_str());

    leg->Draw();


}

int CombPlots(std::string nome,TH3F *map,TH2F *map2)
{
  Int_t steps = 20;
  Int_t size,sizebase;
  size = nome.length();
  sizebase = size - 12;
  char base[120];
  strncpy(base,nome.c_str(),sizebase);
  base[sizebase]='\0';
  //  std::string outfileStr = Form("%s%s", base.c_str(), "ZIP.root");
  char strseed[5];
  strncpy(strseed,(nome.c_str()+sizebase),4);
  strseed[4] = '\0';
  cout << "seed: " << strseed << endl;
  cout << base << endl;
  Int_t seed = atoi(strseed);
  cout << seed << endl;
  std::string buff;
  for (Int_t k = 1 ; k < steps +1 ; k++)
    {
      std::string inFile = Form("%s%dZIP.root",base,seed+k-1);
      cout << "--------------------------------------------------------------------------\nReading: " << inFile.c_str()<<"\n--------------------------------------------------------------------------" << endl;

      TFile* file = TFile::Open(inFile.c_str()); 
      file->cd();
      TNtuple *ntp = (TNtuple*) file->Get("ntp1");
      Int_t n=ntp->GetEntries();
      Float_t phi, rl, zeta;
      ntp->SetBranchAddress("phi",&phi);
      ntp->SetBranchAddress("radlengths",&rl);
      ntp->SetBranchAddress("zeta",&zeta);
      for (Int_t j = 0 ; j < n ; j++)
	{
      	  ntp->GetEntry(j);
	  map->Fill(phi,zeta,rl);  
	  map2->Fill(zeta,rl);  
	}// end for event
    } 
  
  
 

  return 0;
}