diff --git a/.vscode/settings.json b/.vscode/settings.json index 4798fc2..b8d11ef 100644 --- a/.vscode/settings.json +++ b/.vscode/settings.json @@ -1,13 +1,5 @@ { "files.associations": { - "AutoFit.C": "cpp", - "peaks.C": "cpp", - "array": "cpp", - "deque": "cpp", - "list": "cpp", - "vector": "cpp", - "string_view": "cpp", - "span": "cpp", - "SaveTH1IntoText.C": "cpp" + "*.C": "cpp" } } \ No newline at end of file diff --git a/caliPeaks.cpp b/caliPeaks.cpp index 5c78fd5..9e03bf1 100644 --- a/caliPeaks.cpp +++ b/caliPeaks.cpp @@ -65,8 +65,8 @@ void caliPeaks(){ // } // } - TFile * file = new TFile("R9-172_summed_tree.root"); - + TString inFileName = "R9-172_summed_tree.root"; + TFile * file = new TFile(inFileName); TTree * tree = (TTree *) file->Get("rxtree"); const int nBin = (ExRange[1] - ExRange[0])/ExRel; @@ -106,8 +106,8 @@ void caliPeaks(){ Double_t thetaCM; if( SAVETREE ){ - TFile * newFile = new TFile("Cali_132Sn.root", "recreate"); - TTree * newTree = new TTree ("tree", "Calibrated Ex"); + newFile = new TFile("Cali_132Sn_pal.root", "recreate"); + newTree = new TTree ("tree", "Calibrated Ex"); newTree->Branch("Ex_org", &Ex_org, "Ex_org/D"); newTree->Branch("mod", &mod_C, "mod/I"); @@ -122,6 +122,8 @@ void caliPeaks(){ const double ya = 0.854; const double y2 = 2.005; + tree->SetImplicitMT(true); + while( reader.Next() ){ if( *ebis == 0 ) continue; @@ -131,7 +133,7 @@ void caliPeaks(){ int r = *row -1; int p = *pStrip; - // printf("%d %d %d | %f\n", m, r, p, *Ex/1000.); + if( reader.GetCurrentEntry() < 100) printf("%llu | %d %d %d | %f\n", reader.GetCurrentEntry(), m, r, p, *Ex/1000.); if( peaks[m][r][p].size() == 0 ) continue; diff --git a/parallel.cpp b/parallel.cpp new file mode 100644 index 0000000..e1ea79a --- /dev/null +++ b/parallel.cpp @@ -0,0 +1,233 @@ +#include "TFile.h" +#include "TTree.h" +#include "TH1F.h" +#include "TStyle.h" +#include "TH2F.h" +#include "TCanvas.h" +#include "TTreeReader.h" + +#include "AutoFit.C" +#include "SaveTH1IntoText.C" + +#include "TROOT.h" +#include "ROOT/TThreadedObject.hxx" +#include "ROOT/TTreeProcessorMT.hxx" + +// mod = 0, 1, 2 +// row = 1, 2, 3 +// pstrip = 0 - 127 + +#define MAXMOD 3 +#define MAXROW 3 +#define MAXPStrip 128 + +const double ExRange[2] = {-0.5, 3.5}; +const double ExRel = 0.05; //MeV +const int nBin = (ExRange[1] - ExRange[0])/ExRel; + +#define SAVETREE false + +//^==================================== +TH1F * haha[MAXMOD][MAXROW][MAXPStrip]; + +TH2F * dudu[MAXMOD][MAXROW]; + +TH1F * hEx; + +//^==================================== +void parallel(){ + + + /// Load peaks + std::vector peaks[MAXMOD][MAXROW][MAXPStrip]; + + FILE * inList = fopen("fit_result.txt", "r"); + + char buffer[256]; + + while( fgets(buffer, sizeof(buffer), inList) ){ + + std::vector haha = SplitStrAF(buffer, " "); + + // printf("%zu| %s", haha.size(), buffer); + int m = atoi(haha[0].c_str()); + int r = atoi(haha[1].c_str()) - 1; + int p = atoi(haha[2].c_str()); + + for( int i = 3; i < (int) haha.size(); i++ ) peaks[m][r][p].push_back(atof(haha[i].c_str())); + } + + fclose(inList); + + //^=========================== histogram + for( int i = 0; i < MAXMOD; i++ ){ + for( int j = 0; j < MAXROW; j++ ){ + + dudu[i][j] = new TH2F(Form("dudu_%d%d", i, j+1), Form("%d-%d; pStrip; Ex", i, j+1), 128, 0, 128, nBin, ExRange[0], ExRange[1]); + + for( int p = 0; p < MAXPStrip; p ++ ){ + haha[i][j][p] = new TH1F(Form("h_%d%d_%d", i, j+1, p), Form("%d-%d-%d", i, j+1, p), nBin, ExRange[0], ExRange[1]); + } + } + } + + hEx = new TH1F("hEx", "Ex", 10*nBin, ExRange[0], ExRange[1]); + + //^=========================== tree + int nThreads = 4; + ROOT::EnableImplicitMT(nThreads); + + + TString inFileName = "R9-172_summed_tree.root"; + + ROOT::TTreeProcessorMT tp(inFileName, "rxtree", nThreads); + + const double y1 = 0.000; + const double ya = 0.854; + const double y2 = 2.005; + + TFile * newFile = nullptr; + TTree * newTree = nullptr; + + Double_t Ex_org; + Double_t Ex_C; + Int_t mod_C; + Int_t row_C; + Int_t pStrip_C; + Double_t z; + Double_t thetaCM; + + if( SAVETREE ){ + newFile = new TFile("Cali_132Sn_pal.root", "recreate"); + + newTree = new TTree("tree", "Calibred Ex"); + newTree->Branch("Ex_org", &Ex_org, "Ex_org/D"); + newTree->Branch("mod", &mod_C, "mod/I"); + newTree->Branch("row", &row_C, "row/I"); + newTree->Branch("pstrip", &pStrip_C, "pstrip/I"); + newTree->Branch("Ex", &Ex_C, "Ex/D"); + newTree->Branch("z", &z, "z/D"); + newTree->Branch("thetaCM", &thetaCM, "thetaCM/D"); + } + + std::mutex treeMuTex; + + //^========================= define the process in each thread + auto myFunction = [&](TTreeReader &reader) { + + TTreeReaderValue Ex = {reader, "Ex"}; + TTreeReaderValue mod = {reader, "mod"}; + TTreeReaderValue row = {reader, "row"}; + TTreeReaderValue pStrip = {reader, "pstrip"}; + TTreeReaderValue ebis = {reader, "ebis_on"}; + TTreeReaderValue ThetaCM = {reader, "ThetaCM"}; + TTreeReaderValue Zmeasured = {reader, "Zmeasured"}; + + while( reader.Next() ){ + + if( *ebis == 0 ) continue; + if( *row == 0 ) continue; + + int m = *mod; + int r = *row -1; + int p = *pStrip; + + if( reader.GetCurrentEntry() < 100) printf("%llu | %d %d %d | %f\n", reader.GetCurrentEntry(), m, r, p, *Ex/1000.); + + if( peaks[m][r][p].size() == 0 ) continue; + + double slope = 1; + double offset = 0; + + if( peaks[m][r][p].size() == 1 ) { + if( peaks[m][r][p][0] < 1 ){ + offset = -peaks[m][r][p][0]; + }else{ + offset = y2 - peaks[m][r][p][0]; + } + } + if( peaks[m][r][p].size() >= 2 ){ + double x1 = peaks[m][r][p][0]; + double x2 = peaks[m][r][p].back(); + + if( x1 < 0.2 ) { + if( x2 < 1.8 ) { + slope = (ya-y1)/(x2-x1); + }else{ + slope = (y2-y1)/(x2-x1); + } + offset = y1 - slope * x1; + } + if( x1 >= 0.2 ) { + // printf("%d %d %d| x1 : %f, ya : %f \n", m, r, p, x1, ya); + slope = (y2-ya)/(x2-x1); + offset = ya - slope * x1; + } + } + + Ex_org = *Ex/1000.; + mod_C = m; + row_C = r; + pStrip_C = p; + z = *Zmeasured; + thetaCM = *ThetaCM; + + Ex_C = slope * (*Ex/1000.) + offset; + + + treeMuTex.lock(); + + haha[m][r][p]->Fill(Ex_C); + + dudu[m][r]->Fill(p, Ex_C); + + hEx->Fill( *Ex/1000); + + if( SAVETREE ) { + newFile->cd(); + treeMuTex.lock(); + newTree->Fill(); + } + + treeMuTex.unlock(); + + } + + }; + + //^========================= process data + tp.Process(myFunction); + + if( SAVETREE ) { + newFile->cd(); + newTree->Write(); + printf("============== %llu\n", newTree->GetEntries()); + newFile->Close(); + } + + + // gStyle->SetOptStat(""); + + // TCanvas * canvas = new TCanvas("canvas", "canvas", 2100, 1500); + // canvas->Divide(3,3); + + // canvas->cd(1); dudu[0][0]->Draw("colz"); + // canvas->cd(2); dudu[0][1]->Draw("colz"); + // canvas->cd(3); dudu[0][2]->Draw("colz"); + + // canvas->cd(4); dudu[1][0]->Draw("colz"); + // canvas->cd(5); dudu[1][1]->Draw("colz"); + // canvas->cd(6); dudu[1][2]->Draw("colz"); + + // canvas->cd(7); dudu[2][0]->Draw("colz"); + // canvas->cd(8); dudu[2][1]->Draw("colz"); + // canvas->cd(9); dudu[2][2]->Draw("colz"); + + TCanvas * canvasEx = new TCanvas("canvasEx", "canvasEx", 800, 400); + canvasEx->cd(); + hEx->Draw(); + + // TCanvas * chaha = new TCanvas("chaha", "c haah", 500, 500); + // haha[0][1][100]->Draw(); + +} \ No newline at end of file