Skip to content

N
NeuRad_tests
Commit 85a2dca8 authored by Muzalevsky I.A's avatar Muzalevsky I.A
Browse files
Options

script for analizing data with CFD method modified

parent def171cc
Hide whitespace changes
Inline Side-by-side
Showing with 31 additions and 22 deletions
convertTektronix/IntegralCFD.c
View file @ 85a2dca8
void IntegralCFD(){ void IntegralCFD(){
gSystem->Load("../libData.so");
const Int_t Tdelta=12; const Int_t Tdelta=12;
const Double_t coeff=0.5; const Double_t coeff=0.5;
Int_t j,i,nentry,NumM,imax,imin; Int_t j,i,nentry,NumM,imax,imin;
Double_t A[1000]; /// амплитуды с сырых данных (отрицательные) Double_t A1[1000],A2[1000],A3[1000],A4[1000]; /// амплитуды с сырых данных (отрицательные)
Double_t T[1000]; Double_t T[1000];
Double_t Am[1000]; // амплитуды с перевёрнутых пиков Double_t Am[1000]; // амплитуды с перевёрнутых пиков
Double_t Tim[1000]; Double_t Tim[1000];
...@@ -14,23 +15,32 @@ void IntegralCFD(){ ...@@ -14,23 +15,32 @@ void IntegralCFD(){
//TF1 *func = new TF1("fit",p,0,20e-09,1); // функция фитирования //TF1 *func = new TF1("fit",p,0,20e-09,1); // функция фитирования
TF1 *fit1 = new TF1("fit1","[0]"); TF1 *fit1 = new TF1("fit1","[0]");
fit1->SetParName(0,"Rlevel"); fit1->SetParName(0,"Rlevel");
fit1->SetRange(0,2e-08); fit1->SetRange(0,20);
TF1 *fit2 = new TF1("fit2","[0]"); ////////// создание функций фитирования TF1 *fit2 = new TF1("fit2","[0]"); ////////// создание функций фитирования
fit2->SetRange(8e-08,1e-07); fit2->SetRange(80,100);
fit2->SetParName(0,"Llevel"); fit2->SetParName(0,"Llevel");
TF1 *fit3 = new TF1("fit3","[0]*x+[1]"); TF1 *fit3 = new TF1("fit3","[0]*x+[1]");
fit3->SetParName(0,"LinPar"); fit3->SetParName(0,"LinPar");
TFile *f = new TFile("exp2.root","READ"); TFile *f = new TFile("exp2.root","READ");
TTree *tree1 = (TTree*)f->Get("tree"); TTree *tree1 = (TTree*)f->Get("rtree");
////// //////
TBranch *branch1 = tree1->GetBranch("A"); TBranch *branch1 = tree1->GetBranch("A1");
branch1->SetAddress(A); branch1->SetAddress(A1);
TBranch *branch2 = tree1->GetBranch("A2");
branch2->SetAddress(A2);
TBranch *branch3 = tree1->GetBranch("A3");
branch3->SetAddress(A3);
TBranch *branch4 = tree1->GetBranch("A4");
branch4->SetAddress(A4);
TBranch *branch2 = tree1->GetBranch("T"); TBranch *branch5 = tree1->GetBranch("T");
branch2->SetAddress(T); branch5->SetAddress(T);
////// //////
TFile *f1 = new TFile("exp2Integral.root","RECREATE"); // файл для записи TFile *f1 = new TFile("exp2Integral.root","RECREATE"); // файл для записи
TTree *theTree = new TTree("theTree","peak"); TTree *theTree = new TTree("theTree","peak");
...@@ -40,7 +50,7 @@ void IntegralCFD(){ ...@@ -40,7 +50,7 @@ void IntegralCFD(){
theTree->Branch("zeroLEVEL1",&zeroLEVEL1,"zeroLEVEL1/D"); theTree->Branch("zeroLEVEL1",&zeroLEVEL1,"zeroLEVEL1/D");
theTree->Branch("zCFD",&zCFD,"zCFD/D"); theTree->Branch("zCFD",&zCFD,"zCFD/D");
theTree->Branch("tm",&tm,"tm/D"); //// tm - положение пика(максимума амплитуды) записывается время(T[i]) при котором этот максимум был theTree->Branch("tm",&tm,"tm/D"); //// tm - положение пика(максимума амплитуды) записывается время(T[i]) при котором этот максимум был
theTree->Branch("t1_15",&t1_15,"t1 _15/D"); /// время когда амплитуда достигла 20 процентов от максимума( смотрим слева направо) theTree->Branch("t1_15",&t1_15,"t1_15/D"); /// время когда амплитуда достигла 20 процентов от максимума( смотрим слева направо)
theTree->Branch("t2_15",&t2_15,"t2_15/D"); theTree->Branch("t2_15",&t2_15,"t2_15/D");
theTree->Branch("m",&m,"m/D"); // максимальная амплитуда theTree->Branch("m",&m,"m/D"); // максимальная амплитуда
theTree->Branch("m1_15",&m1_15,"m1_15/D"); // амплитуда в первой точке превышения 20% theTree->Branch("m1_15",&m1_15,"m1_15/D"); // амплитуда в первой точке превышения 20%
...@@ -64,10 +74,10 @@ void IntegralCFD(){ ...@@ -64,10 +74,10 @@ void IntegralCFD(){
for(i=0;i<1000;i++) {Acfd[i]=0;} for(i=0;i<1000;i++) {Acfd[i]=0;}
for(i=0;i<1000;i++)// нахожу интеграл, полный, точку максимума для фрейма номер j for(i=0;i<1000;i++)// нахожу интеграл, полный, точку максимума для фрейма номер j
{ {
Am[i] = -A[i]; Am[i] = -A1[i];
Tim[i] = T[i]; Tim[i] = T[i];
if( i>Tdelta ){ if( i>Tdelta ){
Acfd[i] = A[i]*coeff; Acfd[i] = A1[i]*coeff;
Acfd[i] = Acfd[i] + Am[i-Tdelta]; Acfd[i] = Acfd[i] + Am[i-Tdelta];
} }
if((Am[i]>m) || (Am[i]==m)){m = Am[i];tm = T[i];} if((Am[i]>m) || (Am[i]==m)){m = Am[i];tm = T[i];}
...@@ -75,8 +85,8 @@ void IntegralCFD(){ ...@@ -75,8 +85,8 @@ void IntegralCFD(){
} }
maxCFD = 0; minCFD = 0; maxCFD = 0; minCFD = 0;
for(i=0;i<1000;i++){ // поиск максимума и минимума у преобразованного сигнала, для поиска пересечения с нулём for(i=0;i<1000;i++){ // поиск максимума и минимума у преобразованного сигнала, для поиска пересечения с нулём
if(Acfd[i]>maxCFD) {maxCFD = Acfd[i];imax=i;t1=Tim[i];} if(Acfd[i]>maxCFD) {maxCFD = Acfd[i];imax=i;}
if(Acfd[i]<minCFD) {minCFD = Acfd[i];imin=i;t2=Tim[i];} if(Acfd[i]<minCFD) {minCFD = Acfd[i];imin=i;}
} }
/*TGraph *gr3 = new TGraph(1000,Tim,Acfd); /*TGraph *gr3 = new TGraph(1000,Tim,Acfd);
...@@ -84,11 +94,11 @@ void IntegralCFD(){ ...@@ -84,11 +94,11 @@ void IntegralCFD(){
c1->Update(); c1->Update();
cout<<endl<<imin<<" "<<imax<<endl; */ cout<<endl<<imin<<" "<<imax<<endl; */
NumM = (tm - 5e-013)*(1e+10); // номер точки (0,1000) с максимальной амплитудой сигнала NumM = (tm - 5e-04)*(10); // номер точки (0,1000) с максимальной амплитудой сигнала
Integrall=0; Integrall=0;
for(i=0;i<1000;i++) // нахожу интеграл, c границами for(i=0;i<1000;i++) // нахожу интеграл, c границами
{ {
if( (Tim[i]>(tm-3e-9)) && (Tim[i]<(tm+17e-9)) ) {Integrall = Integrall + Am[i];} if( (Tim[i]>(tm-3)) && (Tim[i]<(tm+17)) ) {Integrall = Integrall + Am[i];}
} }
//cout<<endl<<"integral= "<<Integrall<<endl; //cout<<endl<<"integral= "<<Integrall<<endl;
////////////////// //////////////////
...@@ -102,7 +112,7 @@ void IntegralCFD(){ ...@@ -102,7 +112,7 @@ void IntegralCFD(){
t2_15=0;m2_15 = 0; mL_15 = 0; t2_15=0;m2_15 = 0; mL_15 = 0;
for(i=NumM;i>0;i--) for(i=NumM;i>0;i--)
{ {
if((Am[i]-zeroLEVEL1)<0.1*(m-zeroLEVEL1)) {t1_15 = Tim[i+1];m1_15 = A[i+1];break;} if((Am[i]-zeroLEVEL1)<0.1*(m-zeroLEVEL1)) {t1_15 = Tim[i+1];m1_15 = Am[i+1];break;}
} }
for(i=NumM;i>0;i--) for(i=NumM;i>0;i--)
{ {
...@@ -114,25 +124,24 @@ void IntegralCFD(){ ...@@ -114,25 +124,24 @@ void IntegralCFD(){
} }
for(i=NumM;i<1000;i++) for(i=NumM;i<1000;i++)
{ {
if((Am[i]-zeroLEVEL1)<0.15*(m-zeroLEVEL1)) {t2_15 = Tim[i-1];m2_15 = A[i-1];break;} if((Am[i]-zeroLEVEL1)<0.15*(m-zeroLEVEL1)) {t2_15 = Tim[i-1];m2_15 = Am[i-1];break;}
} }
for(i=NumM;i<1000;i++) for(i=NumM;i<1000;i++)
{ {
if((Am[i]-zeroLEVEL1)>0.15*(m-zeroLEVEL1)) {tL_15 = Tim[i-1];mL_15 = A[i-1];} if((Am[i]-zeroLEVEL1)>0.15*(m-zeroLEVEL1)) {tL_15 = Tim[i-1];mL_15 = Am[i-1];}
} }
deltaZero = 100; deltaZero = 100;
for(i=imin;i<imax;i++){ // поиск нуля в сигнале CFD for(i=imin;i<imax;i++){ // поиск нуля в сигнале CFD
if(abs(A[i])<deltaZero) {deltaZero = abs(A[i]); Tcfd = T[i];} if(abs(Acfd[i])<deltaZero) {deltaZero = abs(Acfd[i]); Tcfd = T[i];}
} }
fit3->SetRange(t1_15,t1_90); fit3->SetRange(t1_15,t1_90);
gr1->Fit("fit3","R","goff"); gr1->Fit("fit3","R","goff");
c1->Update(); c1->Update();
ChiSquare = fit3->GetChisquare(); ChiSquare = fit3->GetChisquare();
LinPar = fit3->GetParameter(0); LinPar = fit3->GetParameter(0);
theTree->Fill(); theTree->Fill();
} }
////////// //////////
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment