// compile with:
// g++ -o plotchippedsubsh.exe plotchippedsubsh.cpp `root-config --cflags --glibs`
//
#include "TH1F.h"
#include "TProfile.h"
#include "TMath.h"
#include "TF1.h"
#include "TLegend.h"
#include "TCanvas.h"
 //
#include "TROOT.h"
//#include "TStyle.h"
#include "TTree.h"
#include "TFile.h"
#include "TGraphErrors.h"
#include "TStyle.h"
#include "TLine.h"
#include <string>
#include <cmath>
#include <iostream>
#include <fstream>
#include <vector>
#include <algorithm>

#include <sys/types.h>
#include <sys/stat.h>
#include <unistd.h>


//

using namespace std;

int plotchippedsubsh(int mcm = 0,const int iEvent = 1,int chip = 0, int saltronb = 0, int shaper = 0,  const char *dir = "",const char *filename = "", const char *title = "test")
{
  char name[100];
  int sample = 0;
  int peakthr0 = 600;
  int peakthr1 = 300;
  int peaktype = 0;
  int peaktime;
  int peakfound = 0;
  int peakoffset = 0;
  int polarity = 1;
  int nbpeakstype[3];
  int AdcValue;
  int maxSample;
  unsigned int ErrorCode[16] = {0};
  unsigned int WarningCode[16] = {0};
  int item = 0;
  float peak_min_mean[3] = {0};
  float peak_max_mean[3] = {0};
  float peak_max_rms[3] = {0};
  float ped_min_mean[3] = {0};
  float ped_max_mean[3] = {0};
  float ped_max_rms[3] = {0};

  float ratioLow;
  float ratioHigh;
  Char_t iChan; // Channel number
  Char_t iChip; // Channel number
  UShort_t iClLength; //number of samples
  UShort_t iClTime; // Cluster time
  UInt_t iEventnb;
  Int_t iSamples[1024]; // Samples
  TH1F* hSignal[16];
  TH1F* hPeaktype0[16];
  TH1F* hPeaktype1[16];
  TH1F* hPeaktype2[16];
  TH1F* hPed0[16];
  TH1F* hPed1[16];
  TH1F* hPed2[16];
  TH1F* hRatioLow;
  TH1F* hRatioHigh;

  FILE *limitfile;
  sprintf(name, "%s.txt","plotchippedsub");
  limitfile = fopen(name,"r");
  if (limitfile == NULL) {
    return 4;
  }
  while (fgets(name,sizeof(name),limitfile) != NULL) {
    if (name[0] != '#') {
      switch(item) {
      case 0: sscanf(name,"%f %f %f",&peak_min_mean[0],&peak_max_mean[0],&peak_max_rms[0]); break;
      case 1: sscanf(name,"%f %f %f",&peak_min_mean[1],&peak_max_mean[1],&peak_max_rms[1]); break;
      case 2: sscanf(name,"%f %f %f",&peak_min_mean[2],&peak_max_mean[2],&peak_max_rms[2]); break;
      case 3: sscanf(name,"%f %f %f",&ped_min_mean[0],&ped_max_mean[0],&ped_max_rms[0]); break;
      case 4: sscanf(name,"%f %f %f",&ped_min_mean[1],&ped_max_mean[1],&ped_max_rms[1]); break;
      case 5: sscanf(name,"%f %f %f",&ped_min_mean[2],&ped_max_mean[2],&ped_max_rms[2]); break;
      }
      item++;
    }
  }
  fclose(limitfile);
  if (item == 6) {
    printf("%f %f %f\n",peak_min_mean[0],peak_max_mean[0],peak_max_rms[0]);
    printf("%f %f %f\n",peak_min_mean[1],peak_max_mean[1],peak_max_rms[1]);
    printf("%f %f %f\n",peak_min_mean[2],peak_max_mean[2],peak_max_rms[2]);
    printf("%f %f %f\n",ped_min_mean[0],ped_max_mean[0],ped_max_rms[0]);
    printf("%f %f %f\n",ped_min_mean[1],ped_max_mean[1],ped_max_rms[1]);
    printf("%f %f %f\n",ped_min_mean[2],ped_max_mean[2],ped_max_rms[2]);
  }
  else return 5;
  
  TCanvas *c1 = new TCanvas("c1","Single cluster",600,600);

  sprintf(name, "%s/%s.root",dir,filename);
  TFile *inFile = new TFile(name);
  TTree *tree = (TTree*) inFile->Get("tree");
  if(!tree) {
    cerr << "Tree SALTRO events was not found: " << name << "\n";
    return 2;
  }

  tree->SetBranchAddress("iEventnb",&iEventnb);
  tree->SetBranchAddress("iSamples",&iSamples);
  tree->SetBranchAddress("iChan",&iChan);
  tree->SetBranchAddress("iChip",&iChip);
  tree->SetBranchAddress("iClLength",&iClLength);
  tree->SetBranchAddress("iClTime",&iClTime);

  Int_t nentries = tree->GetEntries(); 
  printf("ENTRIES %d\n",nentries);

  for (Int_t j = 0; j < 16; j++) {
    sprintf(name, "hSignal_%d", j);
    hSignal[j] = ((TH1F *)(gROOT->FindObject(name)));
    if (hSignal[j]) delete hSignal[j];
    hSignal[j] = new TH1F(name, "",1024,-0.5,1023.5);
    hSignal[j]->SetMinimum(0);
    hSignal[j]->SetStats(0);
  }

  for (Int_t j = 0; j < 16; j++) {
    sprintf(name, "hPeaktype0_%d", j);
    hPeaktype0[j] = ((TH1F *)(gROOT->FindObject(name)));
    if (hPeaktype0[j]) delete hPeaktype0[j];
    hPeaktype0[j] = new TH1F(name, "",1024,-0.5,1023.5);
    hPeaktype0[j]->SetMinimum(0);
    //    hPeaktype0[j]->SetMaximum(1024);
  }
  for (Int_t j = 0; j < 16; j++) {
    sprintf(name, "hPeaktype1_%d", j);
    hPeaktype1[j] = ((TH1F *)(gROOT->FindObject(name)));
    if (hPeaktype1[j]) delete hPeaktype1[j];
    hPeaktype1[j] = new TH1F(name, "",1024,-0.5,1023.5);
    hPeaktype1[j]->SetMinimum(0);
    //    hPeaktype1[j]->SetMaximum(1024);
  }
  for (Int_t j = 0; j < 16; j++) {
    sprintf(name, "hPeaktype2_%d", j);
    hPeaktype2[j] = ((TH1F *)(gROOT->FindObject(name)));
    if (hPeaktype2[j]) delete hPeaktype2[j];
    hPeaktype2[j] = new TH1F(name, "",1024,-0.5,1023.5);
    hPeaktype2[j]->SetMinimum(0);
    //    hPeaktype2[j]->SetMaximum(1024);
  }
  for (Int_t j = 0; j < 16; j++) {
    sprintf(name, "hPed0_%d", j);
    hPed0[j] = ((TH1F *)(gROOT->FindObject(name)));
    if (hPed0[j]) delete hPed0[j];
    hPed0[j] = new TH1F(name, "",1024,-0.5,1023.5);
    hPed0[j]->SetMinimum(0);
    //    hPeaktype2[j]->SetMaximum(1024);
  }
  for (Int_t j = 0; j < 16; j++) {
    sprintf(name, "hPed1_%d", j);
    hPed1[j] = ((TH1F *)(gROOT->FindObject(name)));
    if (hPed1[j]) delete hPed1[j];
    hPed1[j] = new TH1F(name, "",1024,-0.5,1023.5);
    hPed1[j]->SetMinimum(0);
    //    hPeaktype2[j]->SetMaximum(1024);
  }
  for (Int_t j = 0; j < 16; j++) {
    sprintf(name, "hPed2_%d", j);
    hPed2[j] = ((TH1F *)(gROOT->FindObject(name)));
    if (hPed2[j]) delete hPed2[j];
    hPed2[j] = new TH1F(name, "",1024,-0.5,1023.5);
    hPed2[j]->SetMinimum(0);
    //    hPeaktype2[j]->SetMaximum(1024);
  }
  
  sprintf(name, "hRatioLow");
  hRatioLow = ((TH1F *)(gROOT->FindObject(name)));
  if (hRatioLow) delete hRatioLow;
  hRatioLow = new TH1F(name, "",200,0.1,2.1);
  hRatioLow->SetTitle("RatioLow");
  hRatioLow->GetXaxis()->SetTitle("Entries");
  hRatioLow->GetYaxis()->SetTitle("RatioLow");
  hRatioLow->SetMinimum(0);
  hRatioLow->SetStats(11110);
  sprintf(name, "hRatioHigh");
  hRatioHigh = ((TH1F *)(gROOT->FindObject(name)));
  if (hRatioHigh) delete hRatioHigh;
  hRatioHigh = new TH1F(name, "",200,0.1,2.1);
  hRatioHigh->SetTitle("RatioHigh");
  hRatioHigh->GetXaxis()->SetTitle("Entries");
  hRatioHigh->GetYaxis()->SetTitle("RatioHigh");
  hRatioHigh->SetMinimum(0);
  hRatioHigh->SetStats(11110);

  gStyle->SetOptStat(1110);
  nbpeakstype[0] =  nbpeakstype[1] = nbpeakstype[2] = 0;
  for (Int_t j = 0; j < nentries; j++) {
    tree->GetEntry(j);
    if (iEventnb == iEvent) {
      for(Int_t l = 0; l < iClLength; l++) {
        hSignal[iChan]->AddBinContent(iClTime-l, iSamples[l]);
      }
    }

    if (iChip == chip ) {
      peaktype = 0;
      peakoffset = 253;
      /*
      if (shaper == 0) peakoffset = 253;
      else if (shaper == 1) peakoffset = 253;
      else if (shaper == 3) peakoffset = 252;
      else if (shaper == 7) peakoffset = 251;
      */
      peaktime = iClTime - peakoffset;
      if (iSamples[peaktime+1] > iSamples[peaktime]) {
	peakoffset--;
	peaktime++;
      }
       
      for(Int_t l = 0; l < 12; l++) {
	//	printf("%2d %2d %4d\n",iChan,l,iSamples[iClTime - 253 - 60 * l + sample]);
	//	hPeaktype1[iChan]->Fill(iSamples[iClTime - 253 - 60 * l + sample]);

	int pedsum = 0;
	for (Int_t p = 5; p <15; p++) {
	  pedsum += iSamples[peaktime + p];
	}
	pedsum /= 10;
	if ( l == 0) {
	  peaktype = 0;
	  peakfound = 0;
	  ratioLow = 0.0;
	  ratioHigh = 0.0;
	  if (polarity > 0) {
	    if (shaper == 7) {
	      ratioLow = float(iSamples[peaktime+1]-pedsum)
		+ float(iSamples[peaktime]-pedsum)
		+ float(iSamples[peaktime-1]-pedsum);
	      maxSample= 0;
	      if (iSamples[peaktime+1] > maxSample)
		maxSample = iSamples[peaktime+1];
	      if (iSamples[peaktime] > maxSample) maxSample = iSamples[peaktime];
	      if (iSamples[peaktime-1] > maxSample) maxSample = iSamples[peaktime-1];
	      if (maxSample > peakthr1) {
		peakfound = 1;
		ratioLow = float(maxSample-pedsum) / ratioLow;
		ratioHigh = ratioLow;
	      }
	    }
	    else if (iSamples[peaktime] > peakthr1)  {
	      peakfound = 1;
	      ratioLow = float(iSamples[peaktime+1]-pedsum)/float(iSamples[peaktime]-pedsum);
	      ratioHigh = float(iSamples[peaktime-1]-pedsum)/float(iSamples[peaktime]-pedsum);
	    }
	  }
	  else {
	    if (((pedsum-iSamples[peaktime]) > 0) && (iSamples[peaktime] < peakthr0))  {
	      peakfound = 1;
	      ratioLow = float(pedsum-iSamples[peaktime+1])/float(pedsum-iSamples[peaktime]);
	      ratioHigh = float(pedsum-iSamples[peaktime-1])/float(pedsum-iSamples[peaktime]);
	    }
	  }
	  hRatioLow->Fill(ratioLow);
	  hRatioHigh->Fill(ratioHigh);
	  if (peakfound == 1) {
	    if (shaper == 0) {
	      if ((ratioHigh < 0.9) && (ratioLow > 0.6)) peaktype = 1;
	      else if ((ratioLow <= 0.6) || (ratioHigh >= 0.9))  peaktype = 2;
	    }
	    else if (shaper == 1) {
	      if ((ratioHigh < 0.8) && (ratioLow > 0.4)) peaktype = 1;
	      else if ((ratioLow <= 0.4) || (ratioHigh >= 0.8))  peaktype = 2;
	    }
	    else if (shaper == 3) {
	      if ((ratioHigh < 0.6) && (ratioLow > 0.2)) peaktype = 1;
	      else if ((ratioLow <= 0.2) || (ratioHigh >= 0.6))  peaktype = 2;
	    }
	    else if (shaper == 7) {
	      if (ratioHigh > 0.85) peaktype = 1;
	      else if (ratioHigh <= 0.85) peaktype = 2;
	    }
	  }
	  nbpeakstype[peaktype]++;
	}
	
	if (polarity > 0) AdcValue = iSamples[peaktime] - pedsum;
	else AdcValue = -(iSamples[peaktime] - pedsum);
	
	if (peaktype == 0) {
	  hPeaktype0[iChan]->Fill(AdcValue);
	  hPed0[iChan]->Fill(pedsum);
	}
	else if (peaktype == 1) {
	  hPeaktype1[iChan]->Fill(AdcValue);
	  hPed1[iChan]->Fill(pedsum);
	}
	else if (peaktype == 2) {
	  hPeaktype2[iChan]->Fill(AdcValue);
	  hPed2[iChan]->Fill(pedsum);
	}
	peaktime -= 60;
      } 
    } 
  }

  printf("PEAKTYPE %d %d %d\n",nbpeakstype[0],nbpeakstype[1],nbpeakstype[2]);
  Int_t oldLevel = gErrorIgnoreLevel;
  gErrorIgnoreLevel = kWarning;

  int l1=7*60+peakoffset - 1;
//252;
  int l2=l1+60;
  int l3=l2+60;
  int x1 = l1-20;
  int x2 = l3+20;
  for (Int_t j = 0; j < 16; j++) {
    c1->cd(iChan+1);
    sprintf(name, "Event %d Channel %d",iEvent,j);
    hSignal[j]->SetTitle(name);
    hSignal[j]->GetXaxis()->SetTitle("Sample #");
    hSignal[j]->GetYaxis()->SetTitle("ADC value");
    hSignal[j]->SetAxisRange(0,1024, "Y");
    hSignal[j]->SetAxisRange(x1,x2, "X");
    hSignal[j]->Draw();
    TLine *line1 = new TLine(l1,0,l1,1024);
    line1->SetLineColor(kBlue);
    line1->Draw();
    TLine *line2 = new TLine(l2,0,l2,1024);
    line2->SetLineColor(kRed);
    line2->Draw();
    TLine *line3 = new TLine(l3,0,l3,1024);
    line3->SetLineColor(kBlue);
    line3->Draw();
    if (j == 0) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-event%d-%s.pdf(",dir,saltronb,saltronb,mcm,chip,iEvent,title);
    //    else if ( j == 15) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-event%d-%s.pdf)",dir,saltronb,saltronb,mcm,chip,iEvent,title);
    else  sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-event%d-%s.pdf",dir,saltronb,saltronb,mcm,chip,iEvent,title);
    c1->SaveAs(name,"pdf");
  }
  hRatioLow->Draw();
  c1->SaveAs(name,"pdf");
  hRatioHigh->Draw();
  sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-event%d-%s.pdf)",dir,saltronb,saltronb,mcm,chip,iEvent,title);
  c1->SaveAs(name,"pdf");

  printf("Peaktype0\n");
  for (Int_t j = 0; j < 16; j++) {
      c1->cd(j+1);
      sprintf(name, "Pulstrain peak0 ch%d %s",j ,title);
      hPeaktype0[j]->SetTitle(name);
      hPeaktype0[j]->GetYaxis()->SetTitle("Entries");
      hPeaktype0[j]->GetXaxis()->SetTitle("ADC value");
      hPeaktype0[j]->Draw();
      if (j == 0) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-peaktype0-%s.pdf(",dir,saltronb,saltronb,mcm,chip,title);
      else if ( j == 15) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-peaktype0-%s.pdf)",dir,saltronb,saltronb,mcm,chip,title);
      else  sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-peaktype0-%s.pdf",dir,saltronb,saltronb,mcm,chip,title);
      c1->SaveAs(name,"pdf");
      if ((hPeaktype0[j]->GetMean() < peak_min_mean[0]) || (hPeaktype0[j]->GetMean() > peak_max_mean[0])) ErrorCode[j] |= 0x0001;
      if (hPeaktype0[j]->GetRMS() > peak_max_rms[0]) ErrorCode[j] |= 0x0002;
      printf("%2d\t%6.2f\t%5.3f\t%6.1f\n",j+iChip*16,hPeaktype0[j]->GetMean(),hPeaktype0[j]->GetRMS(),hPeaktype0[j]->GetEntries());
  }
  
  printf("Peaktype1\n");
  for (Int_t j = 0; j < 16; j++) {
      c1->cd(j+1);
      sprintf(name, "Pulstrain peak1 ch%d %s",j ,title);
      hPeaktype1[j]->SetTitle(name);
      hPeaktype1[j]->GetYaxis()->SetTitle("Entries");
      hPeaktype1[j]->GetXaxis()->SetTitle("ADC value");
      hPeaktype1[j]->Draw();
      if (j == 0) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-peaktype1-%s.pdf(",dir,saltronb,saltronb,mcm,chip,title);
      else if ( j == 15) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-peaktype1-%s.pdf)",dir,saltronb,saltronb,mcm,chip,title);
      else  sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-peaktype1-%s.pdf",dir,saltronb,saltronb,mcm,chip,title);
      c1->SaveAs(name,"pdf");
      if ((hPeaktype1[j]->GetMean() < peak_min_mean[1]) || (hPeaktype1[j]->GetMean() > peak_max_mean[1])) {
        ErrorCode[j] |= 0x0004;
      }
      if (hPeaktype1[j]->GetRMS() > peak_max_rms[1]) ErrorCode[j] |= 0x0008;
      printf("%2d\t%6.2f\t%5.3f\t%6.1f\n",j+iChip*16,hPeaktype1[j]->GetMean(),hPeaktype1[j]->GetRMS(),hPeaktype1[j]->GetEntries());
  }

  printf("Peaktype2\n");
  for (Int_t j = 0; j < 16; j++) {
      c1->cd(j+1);
      sprintf(name, "Pulstrain peak2 ch%d %s",j ,title);
      hPeaktype2[j]->SetTitle(name);
      hPeaktype2[j]->GetYaxis()->SetTitle("Entries");
      hPeaktype2[j]->GetXaxis()->SetTitle("ADC value");
      hPeaktype2[j]->Draw();
      if (j == 0) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-peaktype2-%s.pdf(",dir,saltronb,saltronb,mcm,chip,title);
      else if ( j == 15) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-peaktype2-%s.pdf)",dir,saltronb,saltronb,mcm,chip,title);
      else  sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-peaktype2-%s.pdf",dir,saltronb,saltronb,mcm,chip,title);
      c1->SaveAs(name,"pdf");
      if ((hPeaktype2[j]->GetMean() < peak_min_mean[2]) || (hPeaktype2[j]->GetMean() > peak_max_mean[2])) ErrorCode[j] |= 0x0010;
      if (hPeaktype2[j]->GetRMS() > peak_max_rms[2]) ErrorCode[j] |= 0x0020;
      printf("%2d\t%6.2f\t%5.3f\t%6.1f\n",j+iChip*16,hPeaktype2[j]->GetMean(),hPeaktype2[j]->GetRMS(),hPeaktype2[j]->GetEntries());
  }

  printf("Ped0\n");
  for (Int_t j = 0; j < 16; j++) {
      c1->cd(j+1);
      sprintf(name, "Peaktype 0 ped mcm%d saltro%d ch%d %s", mcm,saltronb,j ,title);
      hPed0[j]->SetTitle(name);
      hPed0[j]->GetYaxis()->SetTitle("Entries");
      hPed0[j]->GetXaxis()->SetTitle("ADC value");
      hPed0[j]->Draw();
      if (j == 0) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-ped0-%s.pdf(",dir,saltronb,saltronb,mcm,chip,title);
      else if ( j == 15) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-ped0-%s.pdf)",dir,saltronb,saltronb,mcm,chip,title);
      else  sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-ped0-%s.pdf",dir,saltronb,saltronb,mcm,chip,title);
      c1->SaveAs(name,"pdf");
      if ((hPed0[j]->GetMean() < ped_min_mean[0]) || (hPed0[j]->GetMean() > ped_max_mean[0])) ErrorCode[j] |= 0x0040;
      if (hPed0[j]->GetRMS() > ped_max_rms[0]) ErrorCode[j] |= 0x0080;
      printf("%2d\t%6.2f\t%5.3f\t%6.1f\n",j+iChip*16,hPed0[j]->GetMean(),hPed0[j]->GetRMS(),hPed0[j]->GetEntries());
  }
  printf("Ped1\n");  
  for (Int_t j = 0; j < 16; j++) {
      c1->cd(j+1);
      sprintf(name, "Peaktype 1 ped mcm%d saltro%d ch%d %s", mcm,saltronb,j ,title);
      hPed1[j]->SetTitle(name);
      hPed1[j]->GetYaxis()->SetTitle("Entries");
      hPed1[j]->GetXaxis()->SetTitle("ADC value");
      hPed1[j]->Draw();
      if (j == 0) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-ped1-%s.pdf(",dir,saltronb,saltronb,mcm,chip,title);
      else if ( j == 15) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-ped1-%s.pdf)",dir,saltronb,saltronb,mcm,chip,title);
      else  sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-ped1-%s.pdf",dir,saltronb,saltronb,mcm,chip,title);
      c1->SaveAs(name,"pdf");
      if ((hPed1[j]->GetMean() < ped_min_mean[1]) || (hPed1[j]->GetMean() > ped_max_mean[1])) ErrorCode[j] |= 0x0100;
      if (hPed1[j]->GetRMS() > ped_max_rms[1]) ErrorCode[j] |= 0x0200;
      printf("%2d\t%6.2f\t%5.3f\t%6.1f\n",j+iChip*16,hPed1[j]->GetMean(),hPed1[j]->GetRMS(),hPed1[j]->GetEntries());
  }
  printf("Ped2\n");
  for (Int_t j = 0; j < 16; j++) {
      c1->cd(j+1);
      sprintf(name, "Peaktype 2 ped mcm%d saltro%d ch%d %s", mcm,saltronb,j ,title);
      hPed2[j]->SetTitle(name);
      hPed2[j]->GetYaxis()->SetTitle("Entries");
      hPed2[j]->GetXaxis()->SetTitle("ADC value");
      hPed2[j]->Draw();
      if (j == 0) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-ped2-%s.pdf(",dir,saltronb,saltronb,mcm,chip,title);
      else if ( j == 15) sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-ped2-%s.pdf)",dir,saltronb,saltronb,mcm,chip,title);
      else  sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-ped2-%s.pdf",dir,saltronb,saltronb,mcm,chip,title);
      c1->SaveAs(name,"pdf");
      if ((hPed2[j]->GetMean() < ped_min_mean[2]) || (hPed2[j]->GetMean() > ped_max_mean[2])) ErrorCode[j] |= 0x0400;
      if (hPed2[j]->GetRMS() > ped_max_rms[2]) ErrorCode[j] |= 0x0800;
      printf("%2d\t%6.2f\t%5.3f\t%6.1f\n",j+iChip*16,hPed2[j]->GetMean(),hPed2[j]->GetRMS(),hPed2[j]->GetEntries());
  }
  gErrorIgnoreLevel = oldLevel;
    // Store result in csvfile
  FILE *csvfile;
  sprintf(name,"%s/chip%d/saltro%d-mcm%d-c%d-pedsub-%s.csv",dir,saltronb,saltronb,mcm,chip,title);
  //  sprintf(name, "%s/%s.csv",dir,filename);

  struct stat   fileBuffer;
  int csvexist = 0;
  if (stat (filename, &fileBuffer) == 0) csvexist = 1;
  csvfile = fopen(name,"a");
  if (csvfile == NULL) {
    return 3;
  }

  if (csvexist == 0) {
    fprintf(csvfile,"Channel,Peak0Mean,Peak0RMS,Peak0Entries,");
    fprintf(csvfile,"Peak1Mean,Peak1RMS,Peak1Entries,");
    fprintf(csvfile,"Peak2Mean,Peak2RMS,Peak2Entries,");
    fprintf(csvfile,"Ped0Mean,Ped0RMS,Ped0Entries,");
    fprintf(csvfile,"Ped1Mean,Ped1RMS,Ped1Entries,");
    fprintf(csvfile,"Ped2Mean,Ped2RMS,Ped2Entries,");
    fprintf(csvfile,"Error,Warning\n");   
  }
  
  for (int j = 0; j < 16; j++) {
    fprintf(csvfile,"%2d,%6.2f,%6.3f,%6.1f,",j+iChip*16,hPeaktype0[j]->GetMean(),hPeaktype0[j]->GetRMS(),hPeaktype0[j]->GetEntries());
    fprintf(csvfile,"%6.2f,%6.3f,%6.1f,",hPeaktype1[j]->GetMean(),hPeaktype1[j]->GetRMS(),hPeaktype1[j]->GetEntries());
    fprintf(csvfile,"%6.2f,%6.3f,%6.1f,",hPeaktype2[j]->GetMean(),hPeaktype2[j]->GetRMS(),hPeaktype2[j]->GetEntries());
    fprintf(csvfile,"%6.2f,%6.3f,%6.1f,",hPed0[j]->GetMean(),hPed0[j]->GetRMS(),hPed0[j]->GetEntries());
    fprintf(csvfile,"%6.2f,%6.3f,%6.1f,",hPed1[j]->GetMean(),hPed1[j]->GetRMS(),hPed1[j]->GetEntries());
    fprintf(csvfile,"%6.2f,%6.3f,%6.1f,",hPed2[j]->GetMean(),hPed2[j]->GetRMS(),hPed2[j]->GetEntries());
    fprintf(csvfile,"%08X,%08X\n",ErrorCode[j],WarningCode[j]);
  }
  fclose(csvfile);

  return 0;
}

int main (int argc, char *argv[])
{
  int iret = 0;
  if (argc != 9) {
    printf("not enough arguments: argc = %d\nSyntax should be:\n./plotchippedsubsh.exe eventnb saltro mcm chip shaper dir filename title\n",argc);
    return -1;
  }
  int eventnb = atoi(argv[1]);
  int saltronb = atoi(argv[2]);
  int mcm = atoi(argv[3]);
  int chip = atoi(argv[4]);
  int shaper=atoi(argv[5]);
iret = plotchippedsubsh(mcm,eventnb,chip,saltronb,shaper,argv[6],argv[7],argv[8]);
  return iret;
}