ChargeSpektrum.c added

diff --git a/newCOMBI/ChargeSpektrum.c b/newCOMBI/ChargeSpektrum.c
new file mode 100644 (file)
index 0000000..09d4c16
--- /dev/null
+++ b/newCOMBI/ChargeSpektrum.c
@@ -0,0 +1,201 @@
+#include <iostream>
+#include <fstream>
+#include <vector>
+#include <string>
+#include <stdlib.h>
+
+
+#include "TFile.h"
+#include "TTree.h"
+#include "TNtuple.h"
+#include "TBranch.h"
+#include "TString.h"
+#include "TCanvas.h"
+#include "TH1F.h"
+#include "TLegend.h"
+#include "TF1.h"
+#include "TMath.h"
+
+#define NUMPIXELS 25
+#define RIGHT_BOUNDARY 7000
+#define VETO_THRESHOLD 20
+
+#define FIT_FUNC "gaus"
+#define DRAW_FITS true
+
+#define DATAPATH "/d/salt/data_recovered/maps/DennisDoering/root/Analysedata/"
+struct histogram
+{
+ TH1F* Seed;
+ TH1F* Sum;
+ TH1F* Veto;
+ int zahl;
+};
+
+struct histogram GenerateHisto(Int_t inputRun);
+void FitPerform(TH1F*, Int_t);
+void ChargeSpektrum()
+{
+ Int_t inputRun=34009;
+ //Int_t inputRun=341270;
+ struct histogram arrayHisto;
+ arrayHisto=GenerateHisto(inputRun);
+  inputRun=341270;
+ struct histogram arrayHisto2;
+ arrayHisto2=GenerateHisto(inputRun);
+ TCanvas* c1 = new TCanvas();
+
+ arrayHisto.Seed->Draw("");
+  double xmin = arrayHisto.Seed->GetXaxis()->GetXmin();
+double xmax = arrayHisto.Seed->GetXaxis()->GetXmax();
+double factor = 1640./400.;
+arrayHisto.Seed->GetXaxis()->SetLimits(xmin*factor,xmax*factor);
+ arrayHisto2.Seed->Draw("same");
+ double xmin2 = arrayHisto2.Seed->GetXaxis()->GetXmin();
+double xmax2 = arrayHisto2.Seed->GetXaxis()->GetXmax();
+//double factor = 1640./400.;
+double factor2 = 1640./7000.;
+arrayHisto2.Seed->GetXaxis()->SetLimits(xmin2*factor2,xmax2*factor2);
+ //arrayHisto.Sum->Draw("same");
+ //arrayHisto.Veto->Draw("same");
+ cout<<arrayHisto.zahl<<endl;
+}
+
+struct histogram GenerateHisto(Int_t inputRun)
+{
+
+cout<<inputRun<<endl;
+
+  TH1F*        histNtuple = new TH1F(Form("hist%i",inputRun), "Histogram title", 200, 0, RIGHT_BOUNDARY);
+  TH1F*        histNtupleSum=new TH1F(Form("histSum%i",inputRun), "Histogram title", 200, 0, RIGHT_BOUNDARY);
+  TH1F*        histNtupleVeto=new TH1F(Form("histVeto%i",inputRun), "Histogram title", 200, 0, RIGHT_BOUNDARY);
+         // open the file
+         TString path = TString(DATAPATH) + Form("%i_0.root", inputRun);
+         TFile* f = TFile::Open(path);
+         if (f->IsZombie()) 
+         {
+           // if we cannot open the file, print an error messageForm("hist%i",nHistNtuple) and return immediatly
+           printf("Error: cannot open %s\n", path.Data());
+           exit(0);
+         }
+         
+         // get a pointer to the tree
+         TNtuple* hitNtuple = (TNtuple*) f->Get("hit");
+         TNtuple* hitsNtuple = (TNtuple*) f->Get("hits");
+         
+         // one array element and one branch for each pixel
+         Float_t pixel[NUMPIXELS];
+         TBranch* pixelBranch[NUMPIXELS];
+         for (Int_t cnt=0; cnt<NUMPIXELS; cnt++) 
+         {
+           hitNtuple->SetBranchAddress(Form("p%i",cnt+1), &pixel[cnt], &pixelBranch[cnt]);
+         }
+        UInt_t frame;
+         TBranch* frameBranch;
+         hitNtuple->SetBranchAddress("frame", &frame, &frameBranch);
+        UInt_t seedPixel;
+         TBranch* seedPixelBranch;
+         hitNtuple->SetBranchAddress("pixel", &seedPixel, &seedPixelBranch);
+         
+         // loop over all hits to fill histogram
+         Int_t nentries = hitNtuple->GetEntries();
+         // Double_t meanNotSeedSum=0;
+         // Int_t nhits=0;
+       
+         for (Int_t cnt=0; cnt<nentries; cnt++) 
+         {
+           hitNtuple->GetEntry(cnt);
+           
+           // take only frames with a maximum of one hit
+           // hitsNtuple->GetEntry((Int_t) frame);
+           // if (frameCounts > 1)
+             // continue;
+           
+           
+           // histogram with the single pixel
+           histNtuple->Fill(pixel[12]);
+           
+           // histogram with the summed pixels
+           Double_t pixelSum = 0;
+           for (Int_t i=0; i<NUMPIXELS; i++)
+               {
+             pixelSum += pixel[i];
+           }
+           histNtupleSum->Fill(pixelSum);
+           
+           // histogram with the single pixel and with "veto trigger"
+           // for the "histNtupleVeto" take only hits where only the seed pixel contains charge
+           Double_t notSeedSum = 0;
+           for (Int_t i=0; i<12; i++)
+             notSeedSum += pixel[i];
+           for (Int_t i=13; i<NUMPIXELS; i++)
+             notSeedSum += pixel[i];
+           if (TMath::Abs(notSeedSum) > VETO_THRESHOLD)
+             continue;
+           // meanNotSeedSum+=notSeedSum;
+           // nhits++;
+           histNtupleVeto->Fill(pixel[12]);    // histogram with the single pixel
+           
+         }
+         FitPerform(histNtuple, inputRun);
+         struct histogram arrayHisto;
+         arrayHisto.Seed=histNtuple;
+         arrayHisto.Sum=histNtupleSum;
+         arrayHisto.Veto=histNtupleVeto;
+         arrayHisto.zahl=3;
+return arrayHisto;}
+
+void FitPerform(TH1F* histNtuple, Int_t inputRun)
+{
+histNtuple->Draw("");
+if (DRAW_FITS) 
+       {
+       Float_t posMax = 0;
+       histNtuple->GetXaxis()->SetRange(histNtuple->GetXaxis()->FindBin(50),histNtuple->GetXaxis()->FindBin(RIGHT_BOUNDARY));  // look only for maxima with x greater than 20
+       Int_t xValMax = histNtuple->GetBinCenter(histNtuple->GetMaximumBin());
+        TF1* fitFunc = new TF1("fitFunc",FIT_FUNC,0,RIGHT_BOUNDARY);
+       if (TString(FIT_FUNC)=="gaus") 
+       {
+         histNtuple->Fit(fitFunc, "N,Q,W", "", xValMax-10, xValMax+10);
+         posMax = fitFunc->GetParameter(1);
+         fitFunc->DrawCopy("same");
+       }
+       else if (TString(FIT_FUNC)=="landau") 
+       {
+         Float_t fitMax1 = 1000;
+         Float_t fitMax2 = 1000;
+         Float_t fitMax3 = 1000;
+         Float_t minFitMax = 1000;
+         Float_t maxFitMax = 1000;
+         
+         histNtuple->Fit(fitFunc, "N,Q,W", "", 20, RIGHT_BOUNDARY);
+         fitMax1 = fitFunc->GetParameter(1);
+         fitFunc->DrawCopy("same");
+         histNtuple->Fit(fitFunc, "N,Q,W", "", 20, fitMax1);
+         fitMax2 = fitFunc->GetParameter(1);
+         fitFunc->SetLineColor(kBlue);
+         fitFunc->SetLineStyle(2);     // dashed
+         fitFunc->DrawCopy("same");
+         histNtuple->Fit(fitFunc, "N,Q,W", "", fitMax1, RIGHT_BOUNDARY);
+         fitMax3 = fitFunc->GetParameter(1);
+         fitFunc->SetLineColor(kGreen);
+         fitFunc->DrawCopy("same");
+         fitFunc->SetLineStyle(1);     // normal for the following fits
+         
+         // Sort the three fits and save error estimation
+         minFitMax = TMath::Min(TMath::Min(fitMax1,fitMax2),fitMax3);
+         maxFitMax = TMath::Max(TMath::Max(fitMax1,fitMax2),fitMax3);
+         posMax = fitMax1 + fitMax2 + fitMax3 - minFitMax - maxFitMax;
+         //fitLandauErrorLeft.push_back(posMax - minFitMax);
+         //fitLandauErrorRight.push_back(maxFitMax - posMax);
+         
+               TString  legendEntry = TString(Form("Run %i: %.1f %.1f  %.1f",inputRun,posMax,minFitMax,maxFitMax ));
+          TLegend* leg = new TLegend(0.5,0.5,0.89,0.89);//(0.6,0.7,0.89,0.89);
+//   leg->SetHeader();//"Legend Title");
+  leg->SetFillStyle(0);
+         leg->AddEntry((TObject*) 0, legendEntry, "");
+         leg->SetTextSize(0.02);
+         leg->Draw();
+         }
+       }
+}
\ No newline at end of file