Index: /trunk/MagicSoft/Mars/mtemp/mpisa/macros/pedestalstudies_pisa.C
===================================================================
--- /trunk/MagicSoft/Mars/mtemp/mpisa/macros/pedestalstudies_pisa.C	(revision 4422)
+++ /trunk/MagicSoft/Mars/mtemp/mpisa/macros/pedestalstudies_pisa.C	(revision 4422)
@@ -0,0 +1,685 @@
+/* ======================================================================== *\
+!
+! *
+! * This file is part of MARS, the MAGIC Analysis and Reconstruction
+! * Software. It is distributed to you in the hope that it can be a useful
+! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
+! * It is distributed WITHOUT ANY WARRANTY.
+! *
+! * Permission to use, copy, modify and distribute this software and its
+! * documentation for any purpose is hereby granted without fee,
+! * provided that the above copyright notice appear in all copies and
+! * that both that copyright notice and this permission notice appear
+! * in supporting documentation. It is provided "as is" without express
+! * or implied warranty.
+! *
+!
+!
+!   Author(s): Markus Gaug, 04/2004 <mailto:markus@ifae.es>
+!
+!   Copyright: MAGIC Software Development, 2000-2004
+!
+!
+\* ======================================================================== */
+//////////////////////////////////////////////////////////////////////////////
+//
+// pedestalstudies_pisa.C
+//
+// macro to study the pedestal and pedestalRMS with the number of FADC 
+// slices summed up. 
+//
+/////////////////////////////////////////////////////////////////////////////////
+const TString pedfile = "/data0/Temp/20040422_23203_P_Mrk421_E.root";
+
+void pedestalstudies_pisa(const TString pedname=pedfile, Int_t pixin=0, Int_t pixout=397,TString f2="",TString f3="",TString f4="",TString f5="")
+{
+
+  Int_t loops = 1;
+  Int_t stepsize = 2;
+
+  gStyle->SetOptStat(1111);
+  gStyle->SetOptFit();
+  
+  TArrayF *hmeandiffinn = new TArrayF(loops);
+  TArrayF *hrmsdiffinn  = new TArrayF(loops);
+  TArrayF *hmeandiffout = new TArrayF(loops);
+  TArrayF *hrmsdiffout  = new TArrayF(loops);
+  TArrayF *hmeaninn  = new TArrayF(loops);
+  TArrayF *hmeanout  = new TArrayF(loops);
+  TArrayF *hrmsinn   = new TArrayF(loops);
+  TArrayF *hrmsout   = new TArrayF(loops);
+  TArrayF *hmuinn    = new TArrayF(loops);
+  TArrayF *hmuout    = new TArrayF(loops);
+  TArrayF *hsigmainn = new TArrayF(loops);
+  TArrayF *hsigmaout = new TArrayF(loops);
+
+  TArrayF *hmeandiffinnerr = new TArrayF(loops);
+  TArrayF *hrmsdiffinnerr  = new TArrayF(loops);
+  TArrayF *hmeandiffouterr = new TArrayF(loops);
+  TArrayF *hrmsdiffouterr  = new TArrayF(loops);
+  TArrayF *hmeaninnerr  = new TArrayF(loops);
+  TArrayF *hmeanouterr  = new TArrayF(loops);
+  TArrayF *hrmsinnerr   = new TArrayF(loops);
+  TArrayF *hrmsouterr   = new TArrayF(loops);
+  TArrayF *hmuinnerr    = new TArrayF(loops);
+  TArrayF *hmuouterr    = new TArrayF(loops);
+  TArrayF *hsigmainnerr = new TArrayF(loops);
+  TArrayF *hsigmaouterr = new TArrayF(loops);
+
+
+  MStatusDisplay *display = new MStatusDisplay;
+  display->SetUpdateTime(500);
+  display->Resize(850,700);
+      
+  //
+  // Create a empty Parameter List and an empty Task List
+  // The tasklist is identified in the eventloop by its name
+  //
+  MParList  plist;
+  MTaskList tlist;
+  plist.AddToList(&tlist);
+
+  TString tabname="";
+
+  for (Int_t samples=2; samples < stepsize*loops+1; samples=samples+stepsize)
+    {
+
+      plist.Reset();
+      tlist.Reset();
+      
+      //
+      // Now setup the tasks and tasklist for the pedestals:
+      // ---------------------------------------------------
+      //
+      
+      MReadMarsFile read("Events", pedname);
+      read.DisableAutoScheme();
+     
+	Int_t neventi;
+ 	neventi=1000;
+
+	Float_t samplingRate=500;	//Hz
+	if (f2!=""){
+	read.AddFile(f2);
+	neventi += 1000;
+	}
+
+	if (f3!=""){
+	read.AddFile(f3);
+	neventi += 1000;
+	}
+
+	if (f4!=""){
+	read.AddFile(f4);
+	neventi += 1000;
+	}
+
+	if (f5!=""){
+	read.AddFile(f5);
+	neventi += 1000;
+	}
+
+	const Int_t NEvt = neventi;
+	cout << "E" <<neventi<< " NEvt:" << NEvt << endl;
+	Float_t timevt[NEvt];
+
+      MGeomApply      geomapl;
+      //
+      // Set the extraction range higher:
+      //		
+      //MExtractFixedWindow sigcalc;
+      //MExtractSignal3 sigcalc;
+      	MExtractSlidingWindow sigcalc;
+      sigcalc.SetRange(0,samples-1,0,1);
+
+     
+  
+      MPedCalcPedRun pedcalc;
+      pedcalc.SetRange(0,samples-1,0,0);
+      pedcalc.SetWindowSize((Int_t)sigcalc.GetNumHiGainSamples());
+
+      //MPedCalcFromData pedcalc;
+      //pedcalc.SetfHiGainThreshold(100);
+
+      //
+      // Additionally to calculating the pedestals, 
+      // you can fill histograms and look at them
+      //
+      MFillH fill("MHPedestalCam", "MExtractedSignalCam");
+      tabname="PedCam ";
+      tabname += samples;
+      //fill.SetNameTab(Form("%s%2d","PedCam",samples));
+      fill.SetNameTab(tabname);
+
+      tlist.AddToList(&read);
+      tlist.AddToList(&geomapl);
+      tlist.AddToList(&sigcalc);
+      tlist.AddToList(&pedcalc);
+      tlist.AddToList(&fill);
+      
+      const MGeomCamMagic      geomcam;
+      MPedestalCam       pedcam;
+      MBadPixelsCam      badcam;
+      badcam.AsciiRead("badpixels.dat");  
+      
+      MHPedestalCam      hpedcam;
+      MCalibrationPedCam cpedcam;
+      
+      plist.AddToList(&geomcam);
+      plist.AddToList(&pedcam);
+      plist.AddToList(&hpedcam);
+      plist.AddToList(&cpedcam);
+      plist.AddToList(&badcam);
+      
+      //
+      // Create and setup the eventloop
+      //
+      MRawEvtHeader *evtheader;
+      MRawEvtData *fRawEvt;
+      TH1F *hpedevt = new TH1F("hpedevt","Pedestals Drift",NEvt,0,(Float_t)NEvt/samplingRate);
+      Bool_t rc;
+      MEvtLoop evtloop;
+      
+      evtloop.SetParList(&plist);
+      evtloop.SetDisplay(display);
+      //tlist.SetDisplay(display);
+      //
+      // Execute first analysis
+      //
+      //if (!evtloop.Eventloop())
+      //  return;
+      rc = evtloop.PreProcess();
+      
+      //while (tlist.Process())
+	//{
+	//cout << "bleah" << endl;
+	//}      
+//tlist.SetDisplay(display);
+
+       Int_t fNumHiGainSample=15;
+	Float_t peddata[NEvt];
+	Int_t evtan=0;
+
+      if (rc){
+          while (tlist.Process())
+          {
+            fRawEvt=(MRawEvtData *)plist.FindObject("MRawEvtData");
+            MRawEvtPixelIter pixel(fRawEvt);
+
+            Float_t ped=0;
+	    evtan++;
+
+ 	    while (pixel.Next())
+	     {
+
+ 	       const UInt_t idx    = pixel.GetPixelId();
+
+               // lista di pixel da sommare
+               if (idx >pixout || idx < pixin)
+		  continue;
+
+	        Byte_t *ptr = pixel.GetHiGainSamples();
+                const Byte_t *end = ptr + fNumHiGainSample;
+      
+                UInt_t sum = 0;
+	        UInt_t sqr = -1;
+
+	        do
+	          {
+	            sum += *ptr;
+	            //sqr += *ptr * *ptr;
+	          }
+	        while (++ptr != end);
+      
+                 ped += (Float_t)sum;
+
+                  //cout << " Pixel " << idx << " pede: " << ped << " --- ";
+                }
+ 
+               //rc = evtloop.Process(0);
+	       evtheader = (MRawEvtHeader *)plist.FindObject("MRawEvtHeader");
+
+             ped /= pixout-pixin+1;
+		Int_t evtn=evtheader->GetDAQEvtNumber();
+             evtn=evtan;
+
+//test
+	//ped =  150 + sin(((Float_t)evtn/16.-TMath::Floor(evtn/16))*6.28); 
+	//cout << ped << " ";
+	timeped = evtn/samplingRate;
+	timevt[evtn-1]=timeped;
+     
+        hpedevt->Fill(timeped,ped);
+
+            //	cout <<" evt:" << evtn;
+	     // Save pedestal in array for Fourier analysis
+	     peddata[evtn-1]=ped;
+          }
+      }
+      if (!evtloop.PostProcess()){
+          rc = kFALSE;
+          return;
+      }
+
+      TCanvas &p1 = display->AddTab("ped studies");
+      p1.Divide(1,2);
+      p1.cd(1);
+	hpedevt->GetXaxis()->SetTitle("Time (s)");
+//      hpedevt->Draw();
+
+	TGraph *drift = new TGraph(evtan,timevt,peddata);
+	drift->SetMarkerStyle(21);
+	drift->SetMarkerSize(0.25);
+	drift->GetXaxis()->SetTitle("time (s)");
+	drift->GetYaxis()->SetTitle(" pedestal");
+	drift->Draw("APL");
+
+
+      // Fourier analysis
+	// campionamento a 500Hz (1000 evt) => freq. max a 250Hz
+
+	MFFT pro;
+	TArrayF *pedarray = new TArrayF(NEvt,peddata);
+	TArrayF *fft; // array con risultato FFT 
+	fft = pro.RealFunctionFFT(pedarray);
+	Int_t nsample = (Int_t)(fft->GetSize()); 
+	
+	const Int_t NNEvt = nsample;
+
+	Float_t x[10000];
+	Float_t prova;
+	for (int i=0;i<NNEvt;i++)
+	 x[i]=(Float_t)i/(NNEvt/(samplingRate/2)); 
+
+	Float_t *ffty = fft->GetArray();
+        cout << "Szie:" << fft->GetSize();
+	ffty[0]=0.;
+	TGraph *gr = new TGraph(nsample,x,ffty);
+	p1.cd(2);
+	gr->GetXaxis()->SetTitle("Frequency (Hz)");
+	gr->SetMarkerStyle(20);
+	gr->SetMarkerSize(0.4);
+	gr->Draw("AP");		
+
+
+
+//-------------------------------
+
+
+      // Set status lines
+      display->SetStatusLine1(evtloop.GetName());
+      display->SetStatusLine2(rc ? "Done." : "Error!");
+      // Stop automatic update
+      display->StopUpdate();
+      // Reallow context menus
+      display->SetNoContextMenu(kFALSE);
+      // Reallow user to exit window by File menu
+      display->UnLock();
+      // 
+      // Look at one specific pixel, after all the histogram manipulations:
+      //
+      /*
+      MHGausEvents &hpix = hpedcam.GetAverageHiGainArea(0);
+      hpix.DrawClone("fourierevents");
+      
+      MHGausEvents &lpix = hpedcam.GetAverageHiGainArea(1);
+      lpix.DrawClone("fourierevents");
+
+      hpedcam[170].DrawClone("fourierevents");
+  
+      */
+
+      MHCamera dispped0  (geomcam, "Ped;Pedestal",       "Mean per Slice");
+      MHCamera dispped2  (geomcam, "Ped;PedestalRms",    "RMS per Slice");
+      MHCamera dispped4  (geomcam, "Ped;Mean",           "Fitted Mean per Slice");
+      MHCamera dispped6  (geomcam, "Ped;Sigma",          "Fitted Sigma per Slice");
+      MHCamera dispped9  (geomcam, "Ped;DeltaPedMean",   "Rel. Diff. Mean per Slice (Fit-Calc.)");
+      MHCamera dispped11 (geomcam, "Ped;DeltaRmsSigma",  "Rel. Diff. RMS per Slice (Fit-Calc.)");
+  
+      dispped0.SetCamContent(  pedcam, 0);
+      dispped0.SetCamError(    pedcam, 1);
+      dispped2.SetCamContent(  pedcam, 2);
+      dispped2.SetCamError(    pedcam, 3);
+      
+      dispped4.SetCamContent( hpedcam, 0);
+      dispped4.SetCamError(   hpedcam, 1);
+      dispped6.SetCamContent( hpedcam, 2);
+      dispped6.SetCamError(   hpedcam, 3);
+      dispped9.SetCamContent( hpedcam, 5);
+      dispped9.SetCamError(   hpedcam, 6);
+      dispped11.SetCamContent(hpedcam, 8);
+      dispped11.SetCamError(  hpedcam, 9);
+  
+      dispped0.SetYTitle("Calc. Pedestal per slice [FADC counts]");
+      dispped2.SetYTitle("Calc. Pedestal RMS per slice [FADC counts]");
+      dispped4.SetYTitle("Fitted Mean per slice [FADC counts]");
+      dispped6.SetYTitle("Fitted Sigma per slice [FADC counts]");
+      dispped9.SetYTitle("Rel. Diff. Pedestal per slice Fit-Calc [1]");
+      dispped11.SetYTitle("Rel. Diff. Pedestal RMS per slice Fit-Calc [1]");
+
+
+      // Histogram values
+      tabname="MeanRMS ";
+      tabname += samples;
+      //fill.SetNameTab(tabname);
+      TCanvas &b1 = display->AddTab(tabname);
+      b1.Divide(4,3);
+
+      CamDraw(b1,dispped0,1,4,*hmeaninn,*hmeanout,*hmeaninnerr,*hmeanouterr,samples,stepsize);
+      CamDraw(b1,dispped2,2,4,*hrmsinn,*hrmsout,*hrmsinnerr,*hrmsouterr,samples,stepsize);  
+      CamDraw(b1,dispped4,3,4,*hmuinn,*hmuout,*hmuinnerr,*hmuouterr,samples,stepsize);
+      CamDraw(b1,dispped6,4,4,*hsigmainn,*hsigmaout,*hsigmainnerr,*hsigmaouterr,samples,stepsize); 
+      
+      tabname="MeanRrmSamples ";
+      tabname += samples;
+      display->SaveAsGIF(3*((samples-1)/stepsize)+2,tabname);
+
+      // Differences
+      tabname="RelDiff ";
+      tabname += samples;
+      TCanvas &c4 = display->AddTab(tabname);
+      c4.Divide(2,3);
+      
+      CamDraw(c4,dispped9,1,2,*hmeandiffinn,*hmeandiffout,*hmeandiffinnerr,*hmeandiffouterr,samples,stepsize);
+      CamDraw(c4,dispped11,2,2,*hrmsdiffinn,*hrmsdiffout,*hrmsdiffinnerr,*hrmsdiffouterr,samples,stepsize); 
+
+      tabname="RelDiffSamples ";
+      tabname += samples;
+      display->SaveAsGIF(3*((samples-1)/stepsize)+3,tabname);
+
+   cout << "Loop N. " <<samples << " completed..." <<endl;
+
+    }
+
+  /*
+  TF1 *logg = new TF1("logg","[1]+TMath::Log(x-[0])",1.,30.,2);
+  logg->SetParameters(1.,3.5);
+  logg->SetParLimits(0,-1.,3.);
+  logg->SetParLimits(1,-1.,7.);
+  logg->SetLineColor(kRed);
+  */
+
+  TCanvas *canvas = new TCanvas("PedstudInner","Pedestal Studies Inner Pixels",600,900);
+  canvas->Divide(2,3);
+  canvas->cd(1);
+
+  TGraphErrors *gmeaninn = new TGraphErrors(hmeaninn->GetSize(),
+                                            CreateXaxis(hmeaninn->GetSize(),stepsize),hmeaninn->GetArray(),
+                                            CreateXaxisErr(hmeaninnerr->GetSize(),stepsize),hmeaninnerr->GetArray());
+  gmeaninn->Draw("A*");
+  gmeaninn->SetTitle("Calculated Mean per Slice Inner Pixels");
+  gmeaninn->GetXaxis()->SetTitle("Nr. added FADC slices");
+  gmeaninn->GetYaxis()->SetTitle("Calculated Mean per slice");
+  //  gmeaninn->Fit("pol0");
+  //  gmeaninn->GetFunction("pol0")->SetLineColor(kGreen);
+  //  //  gmeaninn->Fit(logg);
+
+  canvas->cd(2);
+
+  TGraphErrors *gmuinn = new TGraphErrors(hmuinn->GetSize(),
+                                          CreateXaxis(hmuinn->GetSize(),stepsize),hmuinn->GetArray(),
+                                          CreateXaxisErr(hmuinnerr->GetSize(),stepsize),hmuinnerr->GetArray());
+  gmuinn->Draw("A*");
+  gmuinn->SetTitle("Fitted Mean per Slice Inner Pixels");
+  gmuinn->GetXaxis()->SetTitle("Nr. added FADC slices");
+  gmuinn->GetYaxis()->SetTitle("Fitted Mean per Slice");
+  //  gmuinn->Fit("pol0");
+  //  gmuinn->GetFunction("pol0")->SetLineColor(kGreen);
+  //gmuinn->Fit(logg);
+
+
+  canvas->cd(3);
+
+  TGraphErrors *grmsinn = new TGraphErrors(hrmsinn->GetSize(),
+                                           CreateXaxis(hrmsinn->GetSize(),stepsize),hrmsinn->GetArray(),
+                                           CreateXaxisErr(hrmsinnerr->GetSize(),stepsize),hrmsinnerr->GetArray());
+  grmsinn->Draw("A*");
+  grmsinn->SetTitle("Calculated Rms per Slice Inner Pixels");
+  grmsinn->GetXaxis()->SetTitle("Nr. added FADC slices");
+  grmsinn->GetYaxis()->SetTitle("Calculated Rms per Slice");
+  //  //grmsinn->Fit("pol2");
+  //  //grmsinn->GetFunction("pol2")->SetLineColor(kRed);
+  //  grmsinn->Fit(logg);
+
+  canvas->cd(4);
+
+  TGraphErrors *gsigmainn = new TGraphErrors(hsigmainn->GetSize(),
+                                             CreateXaxis(hsigmainn->GetSize(),stepsize),hsigmainn->GetArray(),
+                                             CreateXaxisErr(hsigmainnerr->GetSize(),stepsize),hsigmainnerr->GetArray());
+  gsigmainn->Draw("A*");
+  gsigmainn->SetTitle("Fitted Sigma per Slice Inner Pixels");
+  gsigmainn->GetXaxis()->SetTitle("Nr. added FADC slices");
+  gsigmainn->GetYaxis()->SetTitle("Fitted Sigma per Slice");
+  //  //  gsigmainn->Fit("pol2");
+  //  //  gsigmainn->GetFunction("pol2")->SetLineColor(kRed);
+  //  gsigmainn->Fit(logg);
+
+  canvas->cd(5);
+
+  TGraphErrors *gmeandiffinn = new TGraphErrors(hmeandiffinn->GetSize(),
+                                                CreateXaxis(hmeandiffinn->GetSize(),stepsize),hmeandiffinn->GetArray(),
+                                                CreateXaxisErr(hmeandiffinnerr->GetSize(),stepsize),hmeandiffinnerr->GetArray());
+  gmeandiffinn->Draw("A*"); 
+  gmeandiffinn->SetTitle("Rel. Difference  Mean per Slice Inner Pixels");
+  gmeandiffinn->GetXaxis()->SetTitle("Nr. added FADC slices");
+  gmeandiffinn->GetYaxis()->SetTitle("Rel. Difference Mean per Slice");
+  //  //gmeandiffinn->Fit("pol2");
+  //  //gmeandiffinn->GetFunction("pol2")->SetLineColor(kBlue);
+  //  gmeandiffinn->Fit(logg);
+
+
+  canvas->cd(6);
+
+  TGraphErrors *grmsdiffinn = new TGraphErrors(hrmsdiffinn->GetSize(),
+                                               CreateXaxis(hrmsdiffinn->GetSize(),stepsize),hrmsdiffinn->GetArray(),
+                                               CreateXaxisErr(hrmsdiffinnerr->GetSize(),stepsize),hrmsdiffinnerr->GetArray());
+  grmsdiffinn->Draw("A*");
+  grmsdiffinn->SetTitle("Rel. Difference Sigma per Slice-RMS Inner Pixels");
+  grmsdiffinn->GetXaxis()->SetTitle("Nr. added FADC slices");
+  grmsdiffinn->GetYaxis()->SetTitle("Rel. Difference Sigma per Slice-RMS");
+  //  //grmsdiffinn->Fit("pol2");
+  //  //grmsdiffinn->GetFunction("pol2")->SetLineColor(kBlue);
+  //  grmsdiffinn->Fit(logg);
+
+  canvas->SaveAs("PedestalStudyInner.root");
+  canvas->SaveAs("PedestalStudyInner.ps");
+
+  TCanvas *canvas2 = new TCanvas("PedstudOut","Pedestal Studies Outer Pixels",600,900);
+  canvas2->Divide(2,3);
+  canvas2->cd(1);
+
+  TGraphErrors *gmeanout = new TGraphErrors(hmeanout->GetSize(),
+                                            CreateXaxis(hmeanout->GetSize(),stepsize),hmeanout->GetArray(),
+                                            CreateXaxisErr(hmeanouterr->GetSize(),stepsize),hmeanouterr->GetArray());
+  gmeanout->Draw("A*");
+  gmeanout->SetTitle("Calculated Mean per Slice Outer Pixels");
+  gmeanout->GetXaxis()->SetTitle("Nr. added FADC slices");
+  gmeanout->GetYaxis()->SetTitle("Calculated Mean per Slice");
+  //  gmeanout->Fit("pol0");
+  //  gmeanout->GetFunction("pol0")->SetLineColor(kGreen);
+  //gmeanout->Fit(logg);
+
+  canvas2->cd(2);
+
+  TGraphErrors *gmuout = new TGraphErrors(hmuout->GetSize(),
+                                          CreateXaxis(hmuout->GetSize(),stepsize),hmuout->GetArray(),
+                                          CreateXaxisErr(hmuouterr->GetSize(),stepsize),hmuouterr->GetArray());
+  gmuout->Draw("A*");
+  gmuout->SetTitle("Fitted Mean per Slice Outer Pixels");
+  gmuout->GetXaxis()->SetTitle("Nr. added FADC slices");
+  gmuout->GetYaxis()->SetTitle("Fitted Mean per Slice");
+  //  gmuout->Fit("pol0");
+  //  gmuout->GetFunction("pol0")->SetLineColor(kGreen);
+  //gmuout->Fit(logg);
+
+  canvas2->cd(3);
+
+  TGraphErrors *grmsout = new TGraphErrors(hrmsout->GetSize(),
+                                           CreateXaxis(hrmsout->GetSize(),stepsize),hrmsout->GetArray(),
+                                           CreateXaxisErr(hrmsouterr->GetSize(),stepsize),hrmsouterr->GetArray());
+  grmsout->Draw("A*");
+  grmsout->SetTitle("Calculated Rms per Slice Outer Pixels");
+  grmsout->GetXaxis()->SetTitle("Nr. added FADC slices");
+  grmsout->GetYaxis()->SetTitle("Calculated Rms per Slice");
+  //  //grmsout->Fit("pol2");
+  //  //grmsout->GetFunction("pol2")->SetLineColor(kRed);
+  //  grmsout->Fit(logg);
+
+  canvas2->cd(4);
+
+  TGraphErrors *gsigmaout = new TGraphErrors(hsigmaout->GetSize(),
+                                             CreateXaxis(hsigmaout->GetSize(),stepsize),hsigmaout->GetArray(),
+                                             CreateXaxisErr(hsigmaouterr->GetSize(),stepsize),hsigmaouterr->GetArray());
+  gsigmaout->Draw("A*");
+  gsigmaout->SetTitle("Fitted Sigma per Slice Outer Pixels");
+  gsigmaout->GetXaxis()->SetTitle("Nr. added FADC slices");
+  gsigmaout->GetYaxis()->SetTitle("Fitted Sigma per Slice");
+  //  //gsigmaout->Fit("pol2");
+  //  //gsigmaout->GetFunction("pol2")->SetLineColor(kRed);
+  //  gsigmaout->Fit(logg);
+
+
+  canvas2->cd(5);
+
+  TGraphErrors *gmeandiffout = new TGraphErrors(hmeandiffout->GetSize(),
+                                                CreateXaxis(hmeandiffout->GetSize(),stepsize),hmeandiffout->GetArray(),
+                                                CreateXaxisErr(hmeandiffouterr->GetSize(),stepsize),hmeandiffouterr->GetArray());
+  gmeandiffout->Draw("A*");
+  gmeandiffout->SetTitle("Rel. Difference  Mean per Slice Outer Pixels");
+  gmeandiffout->GetXaxis()->SetTitle("Nr. added FADC slices");
+  gmeandiffout->GetYaxis()->SetTitle("Rel. Difference Mean per Slice");
+  //  //gmeandiffout->Fit("pol2");
+  //w  //gmeandiffout->GetFunction("pol2")->SetLineColor(kBlue);
+  //  gmeandiffout->Fit(logg);
+
+  canvas2->cd(6);
+
+  TGraphErrors *grmsdiffout = new TGraphErrors(hrmsdiffout->GetSize(),
+                                               CreateXaxis(hrmsdiffout->GetSize(),stepsize),hrmsdiffout->GetArray(),
+                                               CreateXaxisErr(hrmsdiffouterr->GetSize(),stepsize),hrmsdiffouterr->GetArray());
+  grmsdiffout->Draw("A*");
+  grmsdiffout->SetTitle("Rel. Difference Sigma per Slice-RMS Outer Pixels");
+  grmsdiffout->GetXaxis()->SetTitle("Nr. added FADC slices");
+  grmsdiffout->GetYaxis()->SetTitle("Rel. Difference Sigma per Slice-RMS");
+  //  //grmsdiffout->Fit("pol2");
+  //  //grmsdiffout->GetFunction("pol2")->SetLineColor(kBlue);
+  //  grmsdiffout->Fit(logg);
+
+
+  canvas2->SaveAs("PedestalStudyOuter.root");
+  canvas2->SaveAs("PedestalStudyOuter.ps");
+
+
+}
+
+
+void CamDraw(TCanvas &c, MHCamera &cam, Int_t i, Int_t j, TArrayF &a1, TArrayF &a2, 
+             TArrayF &a1err, TArrayF &a2err, Int_t samp, Int_t stepsize)
+{
+
+  c.cd(i);
+  MHCamera *obj1=(MHCamera*)cam.DrawCopy("hist");
+  obj1->SetDirectory(NULL);
+  
+  c.cd(i+j);
+  obj1->Draw();
+  ((MHCamera*)obj1)->SetPrettyPalette();
+
+  c.cd(i+2*j);
+  TH1D *obj2 = (TH1D*)obj1->Projection();
+  obj2->SetDirectory(NULL);
+  
+  //      obj2->Sumw2();
+  obj2->Draw();
+  obj2->SetBit(kCanDelete);
+  
+  const Double_t min   = obj2->GetBinCenter(obj2->GetXaxis()->GetFirst());
+  const Double_t max   = obj2->GetBinCenter(obj2->GetXaxis()->GetLast());
+  const Double_t integ = obj2->Integral("width")/2.5066283;
+  const Double_t mean  = obj2->GetMean();
+  const Double_t rms   = obj2->GetRMS();
+  const Double_t width = max-min;
+  
+  if (rms == 0. || width == 0. )
+    return;
+  
+  TArrayI s0(6);
+  s0[0] = 6;
+  s0[1] = 1;
+  s0[2] = 2;
+  s0[3] = 3;
+  s0[4] = 4;
+  s0[5] = 5;
+  
+  TArrayI inner(1);
+  inner[0] = 0;
+  
+  TArrayI outer(1);
+  outer[0] = 1;
+      
+  // Just to get the right (maximum) binning
+  TH1D *half[2];
+  half[0] = obj1->ProjectionS(s0, inner, "Inner");
+  half[1] = obj1->ProjectionS(s0, outer, "Outer");
+
+  half[0]->SetDirectory(NULL);
+  half[1]->SetDirectory(NULL);
+  
+  for (int i=0; i<2; i++)
+    {
+      half[i]->SetLineColor(kRed+i);
+      half[i]->SetDirectory(0);
+      half[i]->SetBit(kCanDelete);
+      half[i]->Draw("same");
+      half[i]->Fit("gaus","Q+");
+
+      if (i==0)
+        {
+          a1[(samp-1)/stepsize] = half[i]->GetFunction("gaus")->GetParameter(1);
+          a1err[(samp-1)/stepsize] = half[i]->GetFunction("gaus")->GetParError(1);
+          if (a1err[(samp-1)/stepsize] > 3.)
+            a1err[(samp-1)/stepsize] = 1.;
+        }
+     if (i==1)
+       {
+         a2[(samp-1)/stepsize] = half[i]->GetFunction("gaus")->GetParameter(1);
+         a2err[(samp-1)/stepsize] = half[i]->GetFunction("gaus")->GetParError(1);
+          if (a2err[(samp-1)/stepsize] > 3.)
+            a2err[(samp-1)/stepsize] = 1.;
+       }
+    }
+  
+  
+}
+
+// -----------------------------------------------------------------------------
+// 
+// Create the x-axis for the event graph
+//
+Float_t *CreateXaxis(Int_t n, Int_t step)
+{
+
+  Float_t *xaxis = new Float_t[n];
+
+  for (Int_t i=0;i<n;i++)
+    xaxis[i] = 2. + step*i;
+
+  return xaxis;
+                 
+}
+
+// -----------------------------------------------------------------------------
+// 
+// Create the x-axis for the event graph
+//
+Float_t *CreateXaxisErr(Int_t n, Int_t step)
+{
+
+  Float_t *xaxis = new Float_t[n];
+
+  for (Int_t i=0;i<n;i++)
+    xaxis[i] = step/2.;
+
+  return xaxis;
+                 
+}