AEvent.cpp 5.34 KB
Newer Older
1 2 3 4 5 6 7 8 9
/*
 * AEvent.cpp
 *
 *  Created on: Dec 28, 2016
 *      Author: daria
 */

#include "AEvent.h"

10
AEvent::AEvent() : fNPoints(1024) {	//fNPoints is number of points in one event, 1024 or 1000
11 12 13 14 15 16

	Init();
	Reset();

}

17 18 19 20 21 22
AEvent::AEvent(const Int_t npoints) : fNPoints(npoints) {

	Init();
	Reset();
}

23 24
AEvent::~AEvent() {
	// TODO Auto-generated destructor stub
25 26
	delete fGraphSignal;
	delete fGraphCFD;
27
	delete fInputEvent;
28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52
}

void AEvent::SetRawDataFile(const char* inprawfile, const char* treename) {

	TString iFileName = inprawfile;
	TFile *fraw = new TFile(iFileName.Data());
	if ( !fraw->IsOpen() ) {
		Error("SetRawDataFile", "File %s was not opened and won't be processed", iFileName.Data());
	}
	TTree *traw = (TTree*)fraw->Get(treename);
	if (!traw) {
		Error("SetRawDataFile", "Tree %s was not found in file %s", treename, iFileName.Data());
	}
}

void AEvent::ProcessEvent() {

	if (fInputEvent == NULL) {
		Warning("AEvent::ProcessEvent", "Input event wasn't set. Function won't be processed.");
		return;
	}

	const Double_t *amp = fInputEvent->GetAmp();
	const Double_t *time = fInputEvent->GetTime();

53
	for(Int_t j = 0; j < fNPoints; j++) {
54 55 56 57
		fAmpPos[j] = amp[j]*(-1.);
		fTime[j] = time[j];
	}

58
	fZeroLevel = FindZeroLevel();
59
	for(Int_t j = 0; j < fNPoints; j++) {
60
		fAmpPos[j] = fAmpPos[j] - fZeroLevel;
61 62
	}

63
	SetMaxAmplitudes();
64

65
	SetGraphs();
66
	FindFrontProperties();
67
	SetCFD();
68
	SetChargeCFD();
69

70 71 72 73 74 75
	return;

}

void AEvent::Reset() {

76
	for (Int_t i = 0; i < fNPoints; i++) {
77 78
		fAmpPos[i] = 0;
		fTime[i] = 0;
79
		fAmpCFD[i] = 0;
80 81
	}

82 83 84
	fEdgeSlope=0.;
	fTime10=0.;
	fTime90=0.;
85 86
	fAmpMax = 0.;
	fTimeAmpMax = 0.;
87 88
	fTimeCFD = 0.;
	fZeroLevel = 0.;
89
	fChargeCFD = 0.;
90 91 92 93 94 95 96 97 98 99 100 101 102
}

void AEvent::SetInputEvent(RawEvent** event) {

	if (event == 0) {
		Warning("AEvent::SetInputEvent", "Input event was set as 0.");
	}
	fInputEvent = *event;

}

void AEvent::Init() {

103 104 105 106
	fAmpPos.Set(fNPoints);
	fTime.Set(fNPoints);
	fAmpCFD.Set(fNPoints);

107 108
	fGraphSignal = new TGraph();
	fGraphCFD = new TGraph();
109 110
	fInputEvent = 0;

111 112 113
	fCFratio = 0.;
	fCFtimeDelay = 0.;

114 115 116
	fNoiseRangeMin = 0.;
	fNoiseRangeMax = 1.;

117 118 119 120
}

void AEvent::SetGraphs() {

121
	fGraphSignal->Set(fNPoints);
122 123

	for (Int_t i=0; i<fNPoints; i++) {
124
		fGraphSignal->SetPoint(i, fTime[i], fAmpPos[i]);
125 126 127 128
	}

	return;
}
129

130
void AEvent::SetCFD() {
131 132

	Double_t level = 100.;	//is necessary to find cfd amplitude value closest to zero
133 134 135 136 137 138

	fGraphCFD->Set(fNPoints);

	//working variables
	Double_t maxCFD = 0., minCFD = 0.;
	Int_t imax = 0, imin = 0;
139 140 141

	for (Int_t i=0; i<fNPoints; i++) {

142 143 144 145 146
		//CFD method
		if(i>fCFtimeDelay) {
			fAmpCFD[i] = fAmpPos[i]*fCFratio*(-1);
			fAmpCFD[i] = fAmpCFD[i] + fAmpPos[i - fCFtimeDelay];
			fGraphCFD->SetPoint(i, fTime[i], fAmpCFD[i]);
147 148
		}

149
		//point for max CFD amplitude
150
		if(fAmpCFD[i] > maxCFD) {
151 152
			maxCFD = fAmpCFD[i];
			imax = i;
153 154
		}

155
		//point for min CFD amplitude
156
		if(fAmpCFD[i] < minCFD) {
157 158
			minCFD = fAmpCFD[i];
			imin = i;
159 160 161 162 163 164
		}

	}

	//finding "zero" of CFD amplitude
	for(Int_t j = imin; j < imax; j++) {   
165 166
		if(abs(fAmpCFD[j]) < level) {
			level = abs(fAmpCFD[j]);
167
			fTimeCFD = fTime[j];
168
		}
169
	}
Muzalevsky I.A's avatar
Muzalevsky I.A committed
170
}
171 172 173
/////// trying to create method
void AEvent::FindFrontProperties() {

174 175
	Int_t NumM = 0.;

176 177 178 179 180
	if (!fGraphSignal) {
		Warning("AEvent::FindZeroLevel", "Graph was not set");
		return;
	}

181 182
	if(fNPoints == 1000) { NumM = fTimeAmpMax*10.; }	//10 = 1000 points / 100 ns 
	if(fNPoints == 1024) { NumM = fTimeAmpMax*5.; }	//5 = 1024 points / 204.8 ns
183 184 185
			
	for(Int_t i=NumM;i>0;i--)
	{
186
		if( fAmpPos[i]<0.1*fAmpMax ) {fTime10 = fTime[i+1];break;}						
187 188 189
	}
	for(Int_t i=NumM;i>0;i--)
	{
190
		if( fAmpPos[i]<0.9*fAmpMax ) {fTime90 = fTime[i+1];break;}						
191
	}
192

193
	TF1 *fit1 = new TF1("fit1","[1]*x+[0]");	//function for one parameter fitting in the range of pmin-pmax
194 195 196
	fit1->SetRange(fTime10,fTime90);
	fGraphSignal->Fit(fit1,"RQN","goff");
	fEdgeSlope = fit1->GetParameter(1);
197 198 199 200

	delete fit1;
}
////////
201
Double_t AEvent::FindZeroLevel() {
202 203 204

	SetGraphs();
	Double_t correction = 0;
205
	TF1 *fit1 = new TF1("fit1","[0]");	//function for one parameter fitting in the range of pmin-pmax
206
	fit1->SetRange(fNoiseRangeMin,fNoiseRangeMax);
207 208 209 210 211

	if (!fGraphSignal) {
		Warning("AEvent::FindZeroLevel", "Graph was not set");
		return 0;
	}
212

213 214
	fGraphSignal->Fit(fit1,"RQN","goff");
	correction = fit1->GetParameter(0);
215

216 217
	delete fit1;
	return correction;
218
}
219

220 221
void AEvent::SetChargeCFD(Int_t tmin, Int_t tmax) {
	
222 223 224 225
	Double_t integral = 0.;					//voltage
	Double_t time_sig = 0;					//approximate signal duration in seconds
	const Double_t res = 50.; 				//resistance 50 Om
	time_sig = (double)(-tmin + tmax)*(1e-9);
226 227

	for(Int_t i = 0; i < fNPoints; i++) {
228
		if( fTime[i] > (fTimeCFD + tmin) && fTime[i] < (fTimeCFD + tmax) ) {
229 230 231
			integral = integral + fAmpPos[i];
		}
	}
232
	fChargeCFD = integral*time_sig/res;
233 234 235 236 237 238 239
//	cout<<fCharge<<endl;
//	printf("\nIntegral is %f , charge is %lf \n", integral, fCharge);

	return;

}

Muzalevsky I.A's avatar
Muzalevsky I.A committed
240
Double_t AEvent::GetfCFD() {
241
		return fTimeCFD;
Muzalevsky I.A's avatar
Muzalevsky I.A committed
242 243 244 245 246
}

Double_t AEvent::GetOnefTime(Int_t i) {
		return fTime[i];
}
247

Muzalevsky I.A's avatar
Muzalevsky I.A committed
248 249
Double_t AEvent::GetOnefAmpPos(Int_t i) {
		return fAmpPos[i];
250
}
251 252

Double_t AEvent::GetT_10() {
253
		return fTime10;
254 255 256
}

Double_t AEvent::GetT_90() {
257
		return fTime90;
258 259 260
}

Double_t AEvent::GetEdgeSlope() {
261 262
		return fEdgeSlope;
}
263

264 265 266 267 268 269 270 271 272 273 274 275 276
void AEvent::SetMaxAmplitudes() {
	Double_t maxAmp = 0.;
	Double_t maxAmpT = 0.;

	maxAmp = fAmpPos[0];
	for(Int_t j=0; j < fNPoints; j++) {
		if(fAmpPos[j] > maxAmp) {
			maxAmp = fAmpPos[j];
			maxAmpT = fTime[j];
		}
	}
	fAmpMax = maxAmp;
	fTimeAmpMax = maxAmpT;
277

278
}