Commit 170451e3 authored by Muzalevsky I.A's avatar Muzalevsky I.A

program for processing raw data from raw root file added

parent b3ddd8eb
Script for process RAW root files.
1) use command "make" or "make convertRawToAnalyzed" in folder Neurad_tests
2) in folder convertRawToAnalyzed use command: "./convertRawToAnalyzed infile ofile Asize"
where is first parameter (TString infile) - path to the root file with rawData
second parameter (TString ofile) - path to the root file with new data
third parameter (Int_t Asize) - size of the of AEvent class object (1000 or 1024)
Example: "./convertRawToAnalyzed ../data/rawDataTektronix/exp2.root ../data/dataTektronix/exp2.root 1000"
void analyse() #include <fstream>
#include "TString.h"
#include "TFile.h"
#include "TError.h"
#include "TTree.h"
#include <stdio.h>
#include <stdlib.h>
#include "../dataClasses/RawEvent.h"
#include "../dataClasses/AEvent.h"
int main(int argc, char* argv[])
{ {
gSystem->Load("../libData.so"); // gSystem->Load("../libData.so");
if (argc != 4) {
// Tell the user how to run the program
std::cerr << "Usage: " << argv[0] << " [list with input files] [outputfile] [number of points in one event (1000 of 1024)]" << std::endl;
/* "Usage messages" are a conventional way of telling the user
* how to run a program if they enter the command incorrectly.
*/
return 1;
}
// Print the user's name:
TString infiles = argv[1];
TString ofile = argv[2];
TString Asize = argv[3];
TFile *f = new TFile("../data/rawDataDSR4/NeuRad_test_07_1.root"); TFile *f = new TFile(infiles.Data());
TTree *tr = (TTree*)f->Get("rtree"); TTree *tr = (TTree*)f->Get("rtree");
const Int_t noBranches = 4; const Int_t noBranches = 4;
...@@ -19,12 +44,12 @@ void analyse() ...@@ -19,12 +44,12 @@ void analyse()
// tr->SetMakeClass(1); // tr->SetMakeClass(1);
TFile *fw = new TFile("../data/dataDSR4/analysis_07_1.root", "RECREATE"); //create .root file with somehow analyzed data TFile *fw = new TFile(ofile.Data(), "RECREATE"); //create .root file with somehow analyzed data
TTree *tw = new TTree("atree", "title of drs4 analysis tree"); //create analysis tree atree in it TTree *tw = new TTree("atree", "title of drs4 analysis tree"); //create analysis tree atree in it
AEvent *wevent[noBranches]; // pointer to the array (of AEvent class) in which analyzed data for each channel will be put AEvent *wevent[noBranches]; // pointer to the array (of AEvent class) in which analyzed data for each channel will be put
for (Int_t j = 0; j<noBranches; j++) { for (Int_t j = 0; j<noBranches; j++) {
wevent[j] = new AEvent(); wevent[j] = new AEvent(atoi(Asize));
bName.Form("Ach%d.", j); bName.Form("Ach%d.", j);
wevent[j]->SetInputEvent(&revent[j]); //takes raw event from RawEvent wevent[j]->SetInputEvent(&revent[j]); //takes raw event from RawEvent
wevent[j]->SetCFratio(cfRatio); wevent[j]->SetCFratio(cfRatio);
...@@ -54,5 +79,5 @@ void analyse() ...@@ -54,5 +79,5 @@ void analyse()
tw->Write(); tw->Write();
fw->Close(); fw->Close();
return; return 1;
} }
#include "TF1.h"
void IntegralFULL(){
Int_t j,i,nentry,NumM;
Double_t A[1000]; /// амплитуды с сырых данных (отрицательные)
Double_t T[1000];
Double_t Am[1000]; // амплитуды с перевёрнутых пиков
Double_t Tim[1000];
Double_t sum,zeroLEVEL1,zeroLEVEL2,tm,t1_15,m,m1_15,t2_15,m2_15,tL_15,mL_15,Integrall,LinPar,ChiSquare,t_90;
TString var;
//TF1 *func = new TF1("fit",p,0,20e-09,1); // функция фитирования
TF1 *fit1 = new TF1("fit1","[0]");
fit1->SetParName(0,"Rlevel");
fit1->SetRange(0,2e-08);
TF1 *fit2 = new TF1("fit2","[0]"); ////////// создание функций фитирования
fit2->SetRange(8e-08,1e-07);
fit2->SetParName(0,"Llevel");
TF1 *fit3 = new TF1("fit3","[0]*x+[1]");
fit3->SetParName(0,"LinPar");
TFile *f = new TFile("exp2.root","READ");
TTree *tree1 = (TTree*)f->Get("tree");
//////
TBranch *branch1 = tree1->GetBranch("A");
branch1->SetAddress(A);
TBranch *branch2 = tree1->GetBranch("T");
branch2->SetAddress(T);
//////
TFile *f1 = new TFile("exp2Integral.root","RECREATE"); // файл для записи
TTree *theTree = new TTree("theTree","peak");
theTree->Branch("Am",Am,"Am[1000]/D");
theTree->Branch("Tim",Tim,"Tim[1000]/D");
theTree->Branch("sum",&sum,"sum/D");
theTree->Branch("zeroLEVEL1",&zeroLEVEL1,"zeroLEVEL1/D");
theTree->Branch("zeroLEVEL2",&zeroLEVEL2,"zeroLEVEL2/D");
theTree->Branch("tm",&tm,"tm/D"); //// tm - положение пика(максимума амплитуды) записывается время(T[i]) при котором этот максимум был
theTree->Branch("t1_15",&t1_15,"t1_15/D"); /// время когда амплитуда достигла 20 процентов от максимума( смотрим слева направо)
theTree->Branch("t2_15",&t2_15,"t2_15/D");
theTree->Branch("m",&m,"m/D"); // максимальная амплитуда
theTree->Branch("m1_15",&m1_15,"m1_15/D"); // амплитуда в первой точке превышения 20%
theTree->Branch("m2_15",&m2_15,"m2_15/D"); // амплитуда во второй точке превышения 20%
theTree->Branch("tL_15",&tL_15,"tL_15/D"); // время последнего превышения 20%
theTree->Branch("mL_15",&mL_15,"mL_15/D"); //амплитуда в последнем превышении 20%
theTree->Branch("Integrall",&Integrall,"Integrall/D"); //интеграл в диапозоне 20 нс
theTree->Branch("LinPar",&LinPar,"LinPar/D"); //параметр фита прямой границы сигнала
theTree->Branch("ChiSquare",&ChiSquare,"ChiSquare/D"); //хи квадрат
///////// заполнение дерева
TCanvas *c1 = new TCanvas("c","zero line");
j=0;
for(j=0;j<tree1->GetEntries();j++)
{
tree1->GetEntry(j);
sum=0; // обнуление интеграла
m=0; // и максимума
for(i=0;i<1000;i++){ // нахожу интеграл, полный, точку максимума для фрейма номер j
Am[i] = -A[i];
Tim[i] = T[i];
if((Am[i]>m) || (Am[i]==m)){m = Am[i];tm = T[i];}
sum = sum + Am[i];
}
NumM = (tm - 5e-013)*(1e+10); // номер точки (0,1000) с максимальной амплитудой сигнала
Integrall=0;
for(i=0;i<1000;i++){ // нахожу интеграл, c границами
if( (Tim[i]>(tm-3e-9)) && (Tim[i]<(tm+17e-9)) ) {Integrall = Integrall + Am[i];}
}
////////////////// делаю только для первых 12 фреймов
{
TGraph *gr1 = new TGraph(1000,Tim,Am); // поиск нуля
gr1->Draw("Al*");
gr1->Fit("fit1","R");
zeroLEVEL1 = fit1->GetParameter(0);
c1->Update();
gr1->Fit("fit2","R");
zeroLEVEL2 = fit2->GetParameter(0);
c1->Update();
m1_15=0;
t1_15=0; //поиск начальной границы поднятия на 20%
t2_15=0;m2_15 = 0; mL_15 = 0;
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;}
}
for(i=NumM;i>0;i--)
{
if((Am[i]-zeroLEVEL1)<0.9*(m-zeroLEVEL1)) {t_90 = Tim[i+1];break;}
}
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;}
}
for(i=NumM;i<1000;i++)
{
if((Am[i]-zeroLEVEL1)>0.15*(m-zeroLEVEL1)) {tL_15 = Tim[i-1];mL_15 = A[i-1];}
}
fit3->SetRange(t1_15,t_90);
gr1->Fit("fit3","R");
c1->Update();
ChiSquare = fit3->GetChisquare();
LinPar = fit3->GetParameter(0);
}
theTree->Fill();
}
//////////
theTree->Write();
f1->Close();
f->Close();
}
Script for convertation of text files acquired by Tektronix 7354 oscilloscope to ROOT file. Script for convertation of text files acquired by Tektronix 7354 oscilloscope to ROOT file.
1) use command "make" or "make convertTektronix" in folder Neurad_tests
2) in folder convertTektronix use command: "./convertTektronix ../data/rawDataTektronix/infiles.dat ../data/rawDataTektronix/exp2.root"
where is first parameter (../data/rawDataTektronix/infiles.dat) - path to the file with rawData from oscilloscope
second parameter (../data/rawDataTektronix/exp2.root) - path to the root file with rawData created by this programm
void analyseRawData()
{
gSystem->Load("../libData.so");
TFile *f = new TFile("../data/rawDataTektronix/exp2.root");
TTree *tr = (TTree*)f->Get("rtree");
const Int_t noBranches = 4;
TString bName;
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);
tr->SetBranchAddress(bName.Data(), &revent[j]); //read the tree tr with raw data and fill array revent with raw data
}
TH1F *hist1 = new TH1F("hist1", "h1 title", 1000, 0, 100);
TH1F *hist2 = new TH1F("hist2", "h2 title", 1000, 0, 100);
TH1F *hist3 = new TH1F("hist3", "h3 title", 1000, 0, 100);
TH1F *hist4 = new TH1F("hist4", "h4 title", 1000, 0, 100);
TF1 *fit1 = new TF1("fit1","-[0]*exp(-x*[1])");
fit1->SetRange(34,50);
fit1->SetParName(1,"tD");
Long64_t nEntries = tr->GetEntries();
//loop over events
for (j = 0; j < nEntries; j++) {
tr->GetEntry(j);
for(Int_t i =0;i<1000;i++){
//hist->SetBinContent(i,revent[0]->GetAmp(i));
hist1->AddBinContent(i,revent[0]->GetAmp(i));
hist2->AddBinContent(i,revent[1]->GetAmp(i));
hist3->AddBinContent(i,revent[2]->GetAmp(i));
hist4->AddBinContent(i,revent[3]->GetAmp(i));
}
}
TCanvas *c1 = new TCanvas("c1","test",10,10,1000,600);
c1->Divide(2,2);
c1->cd(1);
hist1->Draw();
c1->cd(2);
hist2->Draw();
hist2->Fit("fit1","R","");
c1->cd(3);
hist3->Draw();
c1->cd(4);
hist4->Draw();
}
ch1
../data/rawDataTektronix/ch12016.12.07-02.11.54.dat
../data/rawDataTektronix/ch12016.12.07-08.48.05.dat
../data/rawDataTektronix/ch12016.12.07-10.03.01.dat
ch2
...
...
...
ch3
...
...
...
\ No newline at end of file
void read() {
//Double_t N[5]; // массив для чисел ()
//ifstream myfile;
//Double_t string_to_double( const std::string& s );
//const Int_t par1 = 5;
//const Double_t par2 = 5.34;
// ....
// ....
gSystem->Load("../libData.so");
asdasd();
Double_t A[1000];
Double_t T[1000];
Int_t i,j,n,p,Nchannel;
TString line;
Double_t amp,time;
TFile *f1 = new TFile("exp2.root","RECREATE");
TTree *tree = new TTree("tree","signal");
tree->Branch("A",A,"A[1000]/D");
tree->Branch("T",T,"T[1000]/D");
RawEvent *event1 = new RawEvent();
tree->Bronch("rawEvent1", "RawEvent", &event1);
RawEvent *event2 = new RawEvent();
tree->Bronch("rawEvent2", "RawEvent", &event2);
RawEvent *event3 = new RawEvent();
tree->Bronch("rawEvent3", "RawEvent", &event3);
RawEvent *event4 = new RawEvent();
tree->Bronch("rawEvent4", "RawEvent", &event4);
RawEvent **Nevent;
Nevent = new RawEvent *[4];
Nevent[0] = event1;
Nevent[1] = event2;
Nevent[2] = event3;
Nevent[3] = event4;
/////////
ifstream myfile; // открываю файл с названиями файлов с данными
ifstream myfile1;
myfile.open("../data/rawDataTektronix/infiles.dat");
if (myfile.is_open()) {
TString fileName;
while (!myfile.eof()) {p++;
fileName.ReadLine(myfile);
if ( fileName.IsNull() ) break;
if ( fileName.EqualTo("ch1") ) {Nchannel = 1; }
if ( fileName.EqualTo("ch2") ) {Nchannel = 2; }
if ( fileName.EqualTo("ch3") ) {Nchannel = 3; }
if ( fileName.EqualTo("ch4") ) {Nchannel = 4; }
cout<<fileName.Data()<<" "<< p << " nchannel = "<<Nchannel<<endl;
if ( fileName.BeginsWith("ch") || fileName.BeginsWith("//") ) continue; //есть строки ch1,2
myfile1.open( fileName.Data() );
if (myfile1.is_open()) {
//cout<<" I found data! number "<<endl;
amp=time=0;
while (!myfile1.eof()) {
line.ReadLine(myfile1);
n=i/1000; // if ((i==1002)||(i==900)||(i==2002)||(i==8002)) {cout<<i<<" "<<n<<endl;}
time = 5e-013+ (i-1000*n)*1e-010;
sscanf(line.Data(), "%lf", &amp);
A[i-1000*n] = amp;
T[i-1000*n] = time;
Nevent[Nchannel-1]->SetAmp(amp, i-1000*n);
Nevent[Nchannel-1]->SetTime(time, i-1000*n);
if(i-1000*n ==999) {
tree->Fill();
for(j=0;j<1000;j++){
T[j]=0;
A[j]=0;
}
Nevent[Nchannel-1]->Reset();
}
i++;
}
}
else {
cout<<" asdasdasd "<<endl;
Error("read.cpp", "Some error when opening file %s", fileName.Data());
return;
}
myfile1.close();
}
myfile.close();
}
else {Error("read.cpp", "Some error when opening file","infiles.dat");return;}
////////////////
TString name;
name.Form("../data/rawDataTektronix/ch%d2016.12.07-02.11.54.dat",1);
////////////
ifstream myfile1;
myfile1.open(name.Data());
if (!myfile1.is_open()) {
Error("read.c", "Some error when opening file");
return;
}
i=0;
myfile1.close();
tree->Write();
f1->Close();
}
Double_t asdasd() {
cout << "Function called!!!!!!!!!!!" << endl;
}
...@@ -31,8 +31,8 @@ CONVERTRAWTOANALYZED = $(PWD)/convertRawToAnalyzed ...@@ -31,8 +31,8 @@ CONVERTRAWTOANALYZED = $(PWD)/convertRawToAnalyzed
all: libData.so \ all: libData.so \
read_binary_DRS4 \ read_binary_DRS4 \
convertTektronix convertTektronix \
convertRawToAnalyzed
#ROOT html documentation, it will be done as a program which will be alsa compiled by this makefile, program will be as a last condition after all of the libraries #ROOT html documentation, it will be done as a program which will be alsa compiled by this makefile, program will be as a last condition after all of the libraries
#htmldoc: libAculData.so #htmldoc: libAculData.so
...@@ -69,6 +69,14 @@ read_binary_DRS4: $(CONVERTDRS4)/read_binary.cpp libData.so ...@@ -69,6 +69,14 @@ read_binary_DRS4: $(CONVERTDRS4)/read_binary.cpp libData.so
$(CC) -o $(CONVERTDRS4)/read_binary_DRS4 $(CONVERTDRS4)/read_binary.cpp -lm `root-config --cflags --libs` -L $(PWD) -lData -Wl,-rpath,$(PWD) $(CC) -o $(CONVERTDRS4)/read_binary_DRS4 $(CONVERTDRS4)/read_binary.cpp -lm `root-config --cflags --libs` -L $(PWD) -lData -Wl,-rpath,$(PWD)
@echo 'Finished building target: $@' @echo 'Finished building target: $@'
@echo ' ' @echo ' '
convertTektronix: $(CONVERTTEKTRONIX)/read1.cpp libData.so
@echo 'Building target: $@'
@echo 'Invoking: GCC C++ Linker'
# $(CC) -L $(ROOTLIBS) -shared -o"libData.so" $(DATAOBJS) $(DATALIBS)
$(CC) -o $(CONVERTTEKTRONIX)/convertTektronix $(CONVERTTEKTRONIX)/read1.cpp -lm `root-config --cflags --libs` -L $(PWD) -lData -Wl,-rpath,$(PWD)
@echo 'Finished building target: $@'
@echo ' '
convertRawToAnalyzed: $(CONVERTRAWTOANALYZED)/analyse.cpp libData.so convertRawToAnalyzed: $(CONVERTRAWTOANALYZED)/analyse.cpp libData.so
@echo 'Building target: $@' @echo 'Building target: $@'
......
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