integralFormSignal.cpp 3.41 KB
Newer Older
1 2 3
{
	gSystem->Load("../libData.so");

4
	TFile *f = new TFile("../data/dataTektronix/analysisExp7.root");
5 6 7
	TTree *tr = (TTree*)f->Get("atree");

	const Int_t noBranches = 4;
8 9
	Double_t ZeroTime[4],mh[4],RisTime[4],decayT[4];
	Int_t iZero[4],PosZero[4],deltaT[4],Eh[4],Bh[4]; 
10 11 12
	TString bName;
	AEvent *aevent[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++) {
13
		aevent[j] = new AEvent(1000);	//each raw event element is of class RawEvent()
14 15
		bName.Form("Ach%d.", j);
		tr->SetBranchAddress(bName.Data(), &aevent[j]);	//read the tree tr with raw data and fill array revent with raw data
16
		//cout << tr->SetBranchAddress(bName.Data(), &aevent[j]) << endl;	//read the tree tr with raw data and fill array revent with raw data
17
	}
18 19 20 21 22 23 24 25 26 27 28 29 30

	TString hname,Chname;
	TH1F *hist[4];
	for (Int_t i = 0; i < 4; i++) {
		hname.Form("hist%d",i);
		Chname.Form("CHANNEL%d",i);
		hist[i] = new TH1F(hname.Data(), Chname.Data(), 2000, -100, 100);
	}

	TH1F *hist01 = new TH1F("hist11", "h1 raw", 2000, -100, 100);
	TH1F *hist11 = new TH1F("hist21", "h2 raw", 2000, -100, 100);
	TH1F *hist21 = new TH1F("hist31", "h3 raw", 2000, -100, 100);
	TH1F *hist31 = new TH1F("hist41", "h4 raw", 2000, -100, 100);
31 32 33 34 35 36 37

	TF1 *fit1 = new TF1("fit1","-[0]*exp(-x*[1])");
	fit1->SetParName(1,"tD");

	Long64_t nEntries = tr->GetEntries();

	//loop over events
38
	for (Int_t j = 0; j < nEntries; j++) {
39 40 41 42 43 44 45 46 47
		tr->GetEntry(j);
		for(Int_t k=0;k<4;k++) { // loop for channels
			ZeroTime[k] = aevent[k]->GetfCFD(); // get fCFD for channel №k
			iZero[k] = ZeroTime[k]*10; 
			if(j==0) { PosZero[k] = iZero[k]; }
			deltaT[k] = iZero[k] - PosZero[k];
		}
		
			for(Int_t i = 0; i<1000; i++){
48 49 50 51 52 53 54 55 56 57 58
				if( ((i+deltaT[0])>-1) && ((i+deltaT[0])<1000))	hist[0]->AddBinContent(i+500,aevent[0]->GetOnefAmpPos(i+deltaT[0]));	
				hist01->AddBinContent(i+500,aevent[0]->GetOnefAmpPos(i));

				if( ((i+deltaT[1])>-1) && ((i+deltaT[1])<1000))	hist[1]->AddBinContent(i+500,aevent[1]->GetOnefAmpPos(i+deltaT[1]));
				hist11->AddBinContent(i+500,aevent[1]->GetOnefAmpPos(i));

				if( ((i+deltaT[2])>-1) && ((i+deltaT[2])<1000))	hist[2]->AddBinContent(i+500,aevent[2]->GetOnefAmpPos(i+deltaT[2]));	
				hist21->AddBinContent(i+500,aevent[2]->GetOnefAmpPos(i));

				if( ((i+deltaT[3])>-1) && ((i+deltaT[3])<1000))	hist[3]->AddBinContent(i+500,aevent[3]->GetOnefAmpPos(i+deltaT[3]));	
				hist31->AddBinContent(i+500,aevent[3]->GetOnefAmpPos(i));
59 60 61
			}
	}

62 63 64 65 66 67 68 69 70 71
	for(j=0;j<4;j++) {
	//	j=0;
		mh[j] = hist[j]->GetBinContent(hist[j]->GetMaximumBin());

		for(i = hist[j].GetMaximumBin(); i>0;i--) {	// finding the 
			if( hist[j]->GetBinContent(i)<0.9*mh[j] ) {
				Eh[j] = i+1;
				break;
			}
		}
72

73 74 75 76 77 78 79 80
		for(i = hist[j].GetMaximumBin(); i>0;i--) {	// finding the 
			if( hist[j]->GetBinContent(i)<0.1*mh[j] ) {
				Bh[j] = i+1;
				break;
			}
		}
		RisTime[j] = (Eh[j] - Bh[j])*0.1; // время нарастания 
	}
81

82 83 84 85 86 87
	Double_t rMin[4];
	Double_t rMax[4];
	rMin[0]=-10.5; rMax[0] = -5;
	rMin[1]=-15; rMax[1] = -10;
	rMin[2]=-24; rMax[2] = -18;
	rMin[3]=-14.5; rMax[3] = -9;
88 89


90 91 92 93 94 95 96 97 98 99 100 101 102 103
	TCanvas *c1 = new TCanvas("c1","test",10,10,1000,600);
	c1->Divide(2,2);
	for(i=0;i<4;i++) {
		c1->cd(i+1);
		hist[i]->Draw();
		fit1->SetRange(rMin[i],rMax[i]);
		hist[i]->Fit("fit1","R");
		decayT[i] = 1/(fit1->GetParameter(1));
	}
	
	for(i=0;i<4;i++) {
		cout<<decayT[i]<<" that's decay time  of channel number "<<i<<endl;
		cout<<RisTime[i]<< " that's a RiseTime for this channel "<<i<<endl;
	}
104 105

}