Zero level is determined and CFD time is corrected

convertDRS4/read_binary.cpp

@@ -334,7 +334,9 @@ int main(int argc, const char * argv[])
 	rfile->Close();
 
 //	delete event;
+	printf("error1\n");
 	delete[] event;
+	printf("error2\n");
 
 	return 1;
 }

dataClasses/AEvent.cpp

@@ -55,6 +55,14 @@ void AEvent::ProcessEvent() {
 		fTime[j] = time[j];
 	}
 
+//	SetGraphs();
+	const Double_t zeroLevel = FindZeroLevel();
+//	const Double_t zeroLevel = 0;
+	for(Int_t j = 0; j < fNPoints; j++) {
+		//fAmpPos[j] = amp[j]*(-1.) - zeroLevel;
+		fAmpPos[j] = fAmpPos[j] - zeroLevel;
+	}
+
 	Double_t maxAmp = 0.;
 	Double_t maxAmpT = 0.;
 
@@ -125,9 +133,8 @@ void AEvent::SetGraphs() {
 }
 
 void AEvent::SetCFD() {
-//	Double_t c = 0.5; 		//attenuation coefficient
-//	Int_t td = 5;			//time of delay in points
-	Double_t level = 100.;	//todo what is it?
+
+	Double_t level = 100.;	//is necessary to find cfd amplitude value closest to zero
 
 	fGraphCFD->Set(fNPoints);
 
@@ -139,9 +146,7 @@ void AEvent::SetCFD() {
 
 		//CFD method
 		if(i>fCFtimeDelay) {
-//			fAmpCFD[i] = fAmpPos[i]*cfRatio*(-1);
 			fAmpCFD[i] = fAmpPos[i]*fCFratio*(-1);
-//			fAmpCFD[i] = fAmpCFD[i] + fAmpPos[i - timeDelay];
 			fAmpCFD[i] = fAmpCFD[i] + fAmpPos[i - fCFtimeDelay];
 			fGraphCFD->SetPoint(i, fTime[i], fAmpCFD[i]);
 		}
@@ -167,8 +172,30 @@ void AEvent::SetCFD() {
 			fCFD = fTime[j];
 		}
 	}
+}
 
-	return;
+Double_t AEvent::FindZeroLevel(Int_t pmin, Int_t pmax) {
+
+//	fGraphZero->Set(fNPoints);
+//	const Double_t *amp = fInputEvent->GetAmp();
+//	const Double_t *time = fInputEvent->GetTime();
+	SetGraphs();
+	Double_t correction = 0;
+	TF1 *fit1 = new TF1("fit1","[0]");
+	fit1->SetRange(pmin,pmax);
+
+//	Warning("AEvent::FindZeroLevel", "Graph was not set");
+	if (!fGraphSignal) {
+		Warning("AEvent::FindZeroLevel", "Graph was not set");
+		return 0;
+	}
+//	fGraphSignal->Print();
+	fGraphSignal->Fit(fit1,"RQN","goff");
+//	fGraphSignal->Fit(fit1,"QR","goff");
+	correction = fit1->GetParameter(0);
+//	printf("zero level %f\n", correction);
+	delete fit1;
+	return correction;
 
 
 }

dataClasses/AEvent.h

@@ -17,6 +17,7 @@
 #include "TString.h"
 #include "TTree.h"
 #include "TFile.h"
+#include "TF1.h"
 //#include "TMath.h"
 
 
@@ -43,6 +44,7 @@ private:
 
 	TGraph *fGraphSignal;
 	TGraph *fGraphCFD;
+	TGraph *fGraphZero;
 
 	RawEvent *fInputEvent;		//!
 
@@ -57,20 +59,33 @@ public:
 
 	void SetRawDataFile(const char* inprawfile, const char* treename);
 	void ProcessEvent();
-//	void Integral()
 	void SetInputEvent(RawEvent** event);
 
 	void SetCFratio(Double_t ratio) { fCFratio = ratio; };	
+	//CFD set attenuation coefficient
+
 	void SetCFtimeDelay(Double_t timeDelay) { fCFtimeDelay = timeDelay; };	
+	//CFD set time delay (in points)
 
 	void Reset();
 	 //Resets arrays to zeros
 
-	TGraph* GetGraph() {
-		//todo
+	TGraph* GetGraphCFD() {
+
 		return fGraphCFD;
-	//	return gSignal;
 	}
+	//draws CFD graphs
+
+	TGraph* GetGraphSignal() {
+
+		return fGraphSignal;
+	}
+	//draws signal shape graphs
+
+	Double_t FindZeroLevel(Int_t pmin = 10, Int_t pmax = 100);
+	//for zero level correction
+	//one parameter fit between pmin and pmax
+	//returns fit parameter i.e. number on which amplitude should be corrected
 
 private:
 	void Init();

macros/analyse.C

@@ -13,7 +13,7 @@ void analyse()
 	RawEvent *revent[noBranches];	// pointer to the array (of RawEvent class) in which raw data for each channel will be put
 	for (Int_t j = 0; j<noBranches; j++) {
 		revent[j] = new RawEvent();	//each raw event element is of class RawEvent()
-		bName.Form("ch%d", j);
+		bName.Form("ch%d.", j);
 		tr->SetBranchAddress(bName.Data(), &revent[j]);	//read the tree tr with raw data and fill array revent with raw data
 	}
 
@@ -37,6 +37,9 @@ void analyse()
 
 	for(Long64_t i = 0; i < nentries; i++) {
 		tr->GetEntry(i);
+		if ( !(tr->GetEntry(i)%100) ) {
+			printf("Found event #%d\n", tr->GetEntry(i));
+		}
 
 		for (Int_t j = 0; j<noBranches; j++) {
 			wevent[j]->Reset();

macros/testShowCFD.cxx

@@ -17,7 +17,7 @@ void testShowCFD()
 	for (Int_t i = 0; i < 10; i++) {
 		gr[i] = 0;
 		tr->GetEntry(i+kFirstEvent);
-		gr[i] = new TGraph(*revent->GetGraph());
+		gr[i] = new TGraph(*revent->GetGraphCFD());
 	}//for over events
 
 

macros/testShowGraphs.cxx

@@ -9,24 +9,33 @@ void testShowGraphs()
 	TTree *tr = (TTree*)fr.Get("atree");
 
 	AEvent *revent = new AEvent();
-	tr->SetBranchAddress("Ach0.",&revent);
+	cout << tr->SetBranchAddress("Ach0.",&revent) << endl;
 
 	TGraph *gr[10];
 
+	cout << "warning" << endl;
+
+	cout << tr->GetEntries() << endl;
+
 	//loop over events
-	for (Int_t i = 0; i < 10; i++) {
+	for (Long64_t i = 0; i < 10; i++) {
 		gr[i] = 0;
-		tr->GetEntry(i+kFirstEvent);
-		gr[i] = new TGraph(*revent->GetGraph());
+	cout << "warning 1" << endl;
+		cout << i+kFirstEvent << endl;
+//		tr->GetEntry(i+kFirstEvent);
+		tr->GetEntry(i);
+	cout << "warning 2" << endl;
+		gr[i] = new TGraph(*revent->GetGraphSignal());
 	}//for over events
 
+	cout << "warning 3" << endl;
 
 	TCanvas *c1 = new TCanvas("c1","test",10,10,1000,600);
 	c1->Divide(3,2);
 
-	for (Int_t i = 0; i < 6; i++) {
-		c1->cd(i+1);
-		gr[i]->Draw();
+	for (Int_t k = 0; k < 6; k++) {
+		c1->cd(k+1);
+		gr[k]->Draw();
 	}
 
 //	c1->cd(2);