From 0895ad57edf55ea6a3e1bd10ab98c7782f320e25 Mon Sep 17 00:00:00 2001 From: "Ryan@SOLARIS_testStation" Date: Mon, 19 Aug 2024 18:53:59 -0400 Subject: [PATCH] SingleHistogram add energy Long for PSD, add SplotPolePlotter_MT.C --- Aux/SplitPolePlotter.C | 7 +- Aux/SplitPolePlotter_MT.C | 194 ++++++++++++++++++++++++++++++++++++++ Aux/script.C | 20 ++-- FSUDAQ.cpp | 5 +- SingleSpectra.cpp | 12 ++- 5 files changed, 225 insertions(+), 13 deletions(-) create mode 100644 Aux/SplitPolePlotter_MT.C diff --git a/Aux/SplitPolePlotter.C b/Aux/SplitPolePlotter.C index f538d1d..387f88a 100644 --- a/Aux/SplitPolePlotter.C +++ b/Aux/SplitPolePlotter.C @@ -1,3 +1,6 @@ +#ifndef SPLITPOLEPLOTTER +#define SPLITPOLEPLOTTER + #include "TFile.h" #include "TChain.h" #include "TH1F.h" @@ -266,4 +269,6 @@ void SplitPolePlotter(TChain *tree, TCutG * pidCut = nullptr, double rhoOffset = canvas->cd(10); hXavg_Theta->Draw("colz"); -} \ No newline at end of file +} + +#endif \ No newline at end of file diff --git a/Aux/SplitPolePlotter_MT.C b/Aux/SplitPolePlotter_MT.C new file mode 100644 index 0000000..5f5dc9f --- /dev/null +++ b/Aux/SplitPolePlotter_MT.C @@ -0,0 +1,194 @@ +#include +#include "TTreeReader.h" +#include "TTreeReaderValue.h" +#include "TTreeReaderArray.h" +#include +#include +#include "ROOT/TProcessExecutor.hxx" +#include "ROOT/TThreadedObject.hxx" + +#include "TH2F.h" +#include "TH1F.h" +#include "TCutG.h" +#include "TCanvas.h" + +#include "SplitPolePlotter.C" +#include "../analyzers/SplitPoleHit.h" + + +void SplotPolePlotter_MT(TChain * chain, const int nThread, TCutG * pidCut = nullptr, double rhoOffset = 0, double rhoScaling = 1, bool isFSUDAQ = true){ + + //^====================== Thread Object, destoryed when merge + ROOT::TThreadedObject pCoin("hCoin", "Coincident ", 16, 0, 16, 16, 0, 16); + ROOT::TThreadedObject phMulti("hMulti", "Multiplicity", 16, 0, 16); + ROOT::TThreadedObject pPID("hPID", "PID; Scin_X ; AnodeB", 200, 0, 20000, 100, 0, isFSUDAQ ? 40000 : 5000); + ROOT::TThreadedObject phXavg_Q("hXavg_Q", "Xavg vs Q ", 200, XMIN, XMAX, 200, 0, isFSUDAQ ? 40000 : 5000); + ROOT::TThreadedObject phF("hF", "Front delay line position", 600, XMIN, XMAX); + ROOT::TThreadedObject phB("hB", "Back delay line position", 600, XMIN, XMAX); + ROOT::TThreadedObject phXavg("hAvg", "Xavg", 600, XMIN, XMAX); + ROOT::TThreadedObject phFocal("hFocal", "Front vs Back ", 200, XMIN, XMAX, 200, XMIN, XMAX); + ROOT::TThreadedObject phXavg_Theta("hXavg_Theta", "Xavg vs Theta ", 200, XMIN, XMAX, 200, 0.5, 1.4); + ROOT::TThreadedObject phRay("hRay", "Ray plot", 400, XMIN, XMAX, 400, -50, 50); + ROOT::TThreadedObject phEx("hEx", "Ex; Ex [MeV]; count/10 keV", 600, -1, 5); + ROOT::TThreadedObject phEx_Multi("hEx_Multi", "Ex vs Multi; Ex; Multi", 600, -1, 5, 16, 0, 16); + + + //^==================== TTreeProcessorMT + ROOT::EnableImplicitMT(nThread); + + std::vector fileList_view; + std::vector fileList; + for( int k = 0; k < chain->GetNtrees(); k++){ + fileList_view.push_back(chain->GetListOfFiles()->At(k)->GetTitle()); + fileList.push_back(chain->GetListOfFiles()->At(k)->GetTitle()); + } + ROOT::TTreeProcessorMT tp(fileList_view, "tree"); + // tp.SetTasksPerWorkerHint(1); + + std::mutex mutex; + + int count = 0; + + //^======================= Define process + auto ProcessTask = [&](TTreeReader &reader){ + + TTreeReaderValue evID = {reader, "evID"}; + TTreeReaderValue multi = {reader, "multi"}; + TTreeReaderArray sn = {reader, "sn"}; + TTreeReaderArray ch = {reader, "ch"}; + TTreeReaderArray e = {reader, "e"}; + TTreeReaderArray e2 = {reader, "e2"}; + TTreeReaderArray e_t = {reader, "e_t"}; + TTreeReaderArray e_f = {reader, "e_f"}; + + mutex.lock(); + count++; + printf("-------------- Thread_ID: %d \n", count); + mutex.unlock(); + + SplitPoleHit hit; + hit.SetMassTablePath("../analyzers/mass20.txt"); + hit.CalConstants("12C", "d", "p", 16, 18); // 80MeV, 5 deg + hit.CalZoffset(0.750); // 1.41 T + + + while (reader.Next()) { + + hit.ClearData(); + + phMulti->Fill(*multi); + + // if( *multi != 9 ) continue; + + for( int i = 0; i < *multi; i++){ + + // t2 = e_t[i]; + // if( t2 < t1 ) printf("entry %lld-%d, timestamp is not in order. %llu, %llu\n", processedEntries, i, t2, t1); + if( i == 0 ) t1 = e_t[i]; + // if( e[i] == 65535 ) continue; + + if( ch[i] == ChMap::ScinR ) {hit.eSR = e[i]; hit.tSR = e_t[i] + e_f[i]/1000;} + if( ch[i] == ChMap::ScinL ) {hit.eSL = e[i]; hit.tSL = e_t[i] + e_f[i]/1000;} + if( ch[i] == ChMap::dFR ) {hit.eFR = e[i]; hit.tFR = e_t[i] + e_f[i]/1000;} + if( ch[i] == ChMap::dFL ) {hit.eFL = e[i]; hit.tFL = e_t[i] + e_f[i]/1000;} + if( ch[i] == ChMap::dBR ) {hit.eBR = e[i]; hit.tBR = e_t[i] + e_f[i]/1000;} + if( ch[i] == ChMap::dBL ) {hit.eBL = e[i]; hit.tBL = e_t[i] + e_f[i]/1000;} + if( ch[i] == ChMap::Cathode ) {hit.eCath = e[i]; hit.tCath = e_t[i] + e_f[i]/1000;} + if( ch[i] == ChMap::AnodeF ) {hit.eAF = e[i]; hit.tAF = e_t[i] + e_f[i]/1000;} + if( ch[i] == ChMap::AnodeB ) {hit.eAB = e[i]; hit.tAB = e_t[i] + e_f[i]/1000;} + + for( int j = i+1; j < sn.GetSize(); j++){ + pCoin->Fill(ch[i], ch[j]); + } + } + + unsigned int dQ = hit.eAB; // delta Q + unsigned int Qt = hit.eSL; // total Q + + if( Qt > 0 && dQ > 0 ) { + pPID->Fill(Qt, dQ); + } + + //=============== PID gate cut + if( pidCut ){ + if( !pidCut->IsInside(Qt, dQ) ) continue; + } + + hit.CalData(2); + + if( hit.theta > 1.2 || 0.5 > hit.theta ) continue; + + if( (!TMath::IsNaN(hit.x1) || !TMath::IsNaN(hit.x2)) ) { + phFocal->Fill(hit.x1, hit.x2); + phF->Fill(hit.x1); + phB->Fill(hit.x2); + phXavg->Fill(hit.xAvg); + phXavg_Q->Fill(hit.xAvg, dQ); + phXavg_Theta->Fill( hit.xAvg, hit.theta); + + for( int i = 0; i < 400; i++){ + double z = -50 + 100/400.*i; + + double x = (z/42.8625 + 0.5)* ( hit.x2-hit.x1) + hit.x1; + + phRay->Fill(x,z); + } + + double ex = hit.Rho2Ex( ((hit.xAvg - rhoOffset)/1000/rhoScaling + hit.GetRho0() ) ); + //if( XMIN < hit.xAvg && hit.xAvg < XMAX) printf("x1 : %6.2f, x2 : %6.2f, xAvg %6.2f cm , ex : %f \n", hit.x1, hit.x2, hit.xAvg, ex); + + phEx->Fill(ex); + phEx_Multi->Fill(ex, *multi); + } + } + + return 0; + }; + + //^============================ Run TP + tp.Process(ProcessTask); + + //^============================ Merge all ThreadedObject + auto coin = pCoin.Merge(); + auto hMulti= phMulti.Merge(); + auto hEx = phEx.Merge(); + auto hEx_Multi = phEx_Multi.Merge(); + auto PID = pPID.Merge(); + auto hFocal = phFocal.Merge(); + auto hF = phF.Merge(); + auto hB = phB.Merge(); + auto hXavg = phXavg.Merge(); + auto hXavg_Q = phXavg_Q.Merge(); + auto hXavg_Theta = phXavg_Theta.Merge(); + auto hRay = phRay.Merge(); + + //^============================== Plot + gStyle->SetOptStat("neiou"); + + TCanvas * canvas = new TCanvas("cc", "Split-Pole", 2500, 1000); + canvas->Divide(5, 2); + + canvas->cd(1); { + PID->DrawCopy("colz"); + if( pidCut ) pidCut->Draw("same"); + } + + canvas->cd(2); hRay->DrawCopy("colz"); + + canvas->cd(3); hF->DrawCopy(); + canvas->cd(4); hB->DrawCopy(); + + canvas->cd(5); hXavg_Q->DrawCopy("colz"); + + canvas->cd(6); hXavg->DrawCopy("colz"); + + canvas->cd(7); hEx->DrawCopy(); + + //canvas->cd(8); coin->DrawCopy("colz"); + canvas->cd(8); hEx_Multi->DrawCopy("colz"); + + canvas->cd(9); canvas->cd(9)->SetLogy(); hMulti->DrawCopy(); + + canvas->cd(10); hXavg_Theta->DrawCopy("colz"); + +} \ No newline at end of file diff --git a/Aux/script.C b/Aux/script.C index af1259d..2cebb08 100644 --- a/Aux/script.C +++ b/Aux/script.C @@ -2,22 +2,24 @@ // #include "../MultiBuilder.cpp" #include "SplitPolePlotter.C" +#include "SplitPolePlotter_MT.C" void script(){ - TChain * chain = new TChain("tree"); + TChain * chain = new TChain("tree"); + // chain->Add("raw_binary/run_13/run013_3000.root"); + // chain->Add("data/run*_3000.root"); + chain->Add("data/12C_dp_*_3000.root"); - chain->Add("raw_binary/run_13/run013_3000.root"); -// chain->Add("data/12C_dp_009_3000.root"); - - TFile * pidCutFile = new TFile("cut_proton.root"); -// TFile * pidCutFile = new TFile("cut_proton_FSU.root"); - TCutG * pidCut = (TCutG *) pidCutFile->Get("protons"); + // TFile * pidCutFile = new TFile("cut_proton.root"); + TFile * pidCutFile = new TFile("cut_proton_FSU.root"); + TCutG * pidCut = (TCutG *) pidCutFile->Get("protons"); - SplitPolePlotter(chain, pidCut, 123.307, 2.75, false); + // SplitPolePlotter(chain, pidCut, 123.307, 2.75, false); // for CoMPASS data + // SplitPolePlotter(chain, pidCut, 123.307, 2.75, true); // faster then MT? + SplotPolePlotter_MT(chain, 5, pidCut, 123.307, 2.75, true); -// SplitPolePlotter(chain, pidCut, 123.307, 2.75, true); //^===================================================== diff --git a/FSUDAQ.cpp b/FSUDAQ.cpp index 0f34b1b..366a697 100644 --- a/FSUDAQ.cpp +++ b/FSUDAQ.cpp @@ -1077,7 +1077,7 @@ void FSUDAQ::UpdateScalar(){ if(digiSettings && digiSettings->isVisible() && digiSettings->GetTabID() == iDigi) digiSettings->UpdateACQStatus(acqStatus); - digiMTX[iDigi].lock(); + // digiMTX[iDigi].lock(); QString blockCountStr = QString::number(digi[iDigi]->GetData()->AggCount); blockCountStr += "/" + QString::number(readDataThread[iDigi]->GetReadCount()); @@ -1105,7 +1105,7 @@ void FSUDAQ::UpdateScalar(){ } } - digiMTX[iDigi].unlock(); + // digiMTX[iDigi].unlock(); } @@ -1810,6 +1810,7 @@ void FSUDAQ::OpenCanvas(){ }else{ canvas->show(); canvas->activateWindow(); + canvas->LoadSetting(); } } diff --git a/SingleSpectra.cpp b/SingleSpectra.cpp index 110bae2..3a4aa8c 100644 --- a/SingleSpectra.cpp +++ b/SingleSpectra.cpp @@ -181,7 +181,11 @@ SingleSpectra::SingleSpectra(Digitizer ** digi, unsigned int nDigi, QString rawD fillHistograms = false; QLabel * lbSettingPath = new QLabel( settingPath , this); - ctrlLayout->addWidget(lbSettingPath, 1, 0, 1, 8); + ctrlLayout->addWidget(lbSettingPath, 1, 0, 1, 6); + + QPushButton * bnSaveButton = new QPushButton("Save Hist. Settings", this); + ctrlLayout->addWidget(bnSaveButton, 1, 6, 1, 2); + connect(bnSaveButton, &QPushButton::click, this, &SingleSpectra::SaveSetting); } @@ -207,6 +211,9 @@ SingleSpectra::SingleSpectra(Digitizer ** digi, unsigned int nDigi, QString rawD for( int j = 0; j < digi[i]->GetNumInputCh(); j++){ if( i < nDigi ) { hist[i][j] = new Histogram1D("Digi-" + QString::number(digi[i]->GetSerialNumber()) +", Ch-" + QString::number(j), "Raw Energy [ch]", nBin, eMin, eMax); + if( digi[i]->GetDPPType() == DPPTypeCode::DPP_PSD_CODE ){ + hist[i][j]->AddDataList("Long Energy", Qt::green); + } }else{ hist[i][j] = nullptr; } @@ -359,6 +366,9 @@ void SingleSpectra::FillHistograms(){ // printf(" ch: %d, last fill idx : %d | %d \n", ch, lastFilledIndex[ID][ch], data); hist[ID][ch]->Fill( data ); + if( digi[i]->GetDPPType() == DPPTypeCode::DPP_PSD_CODE ){ + hist[ID][ch]->Fill( digi[ID]->GetData()->GetEnergy2(ch, lastFilledIndex[ID][ch])); + } hist2D[ID]->Fill(ch, data); } if( histVisibility[ID][ch] ) hist[ID][ch]->UpdatePlot();