diff --git a/NPLib/Detectors/ComptonTelescope/TComptonTelescopeSpectra.cxx b/NPLib/Detectors/ComptonTelescope/TComptonTelescopeSpectra.cxx index 6c362cad0c3513de72858d0c5db6f5d3efbd77ec..061964cc3881342b918ddb7f899cf6f08a851f4b 100644 --- a/NPLib/Detectors/ComptonTelescope/TComptonTelescopeSpectra.cxx +++ b/NPLib/Detectors/ComptonTelescope/TComptonTelescopeSpectra.cxx @@ -86,19 +86,19 @@ void TComptonTelescopeSpectra::InitRawSpectra() // DSSD // FRONT_E_RAW name = "CT"+NPL::itoa(i+1)+"_FRONT_E_RAW"; - AddHisto2D(name, name, fStripX, 1, fStripX+1, 512, 0, 8192, "COMPTONTELESCOPE/RAW/FRONTE"); + AddHisto2D(name, name, fStripX, 1, fStripX+1, 512, 0, 1024, "COMPTONTELESCOPE/RAW/FRONTE"); // BACK_E_RAW name = "CT"+NPL::itoa(i+1)+"_BACK_E_RAW"; - AddHisto2D(name, name, fStripY, 1, fStripY+1, 512, 0, 8192, "COMPTONTELESCOPE/RAW/BACKE"); + AddHisto2D(name, name, fStripY, 1, fStripY+1, 512, 0, 1024, "COMPTONTELESCOPE/RAW/BACKE"); // FRONT_T_RAW name = "CT"+NPL::itoa(i+1)+"_FRONT_T_RAW"; - AddHisto2D(name, name, fStripX, 1, fStripX+1, 512, 0, 8192, "COMPTONTELESCOPE/RAW/FRONTT"); + AddHisto2D(name, name, fStripX, 1, fStripX+1, 512, 0, 1024, "COMPTONTELESCOPE/RAW/FRONTT"); // BACK_T_RAW name = "CT"+NPL::itoa(i+1)+"_BACK_T_RAW"; - AddHisto2D(name, name, fStripY, 1, fStripY+1, 512, 0, 8192, "COMPTONTELESCOPE/RAW/BACKT"); + AddHisto2D(name, name, fStripY, 1, fStripY+1, 512, 0, 1024, "COMPTONTELESCOPE/RAW/BACKT"); // FRONT_RAW_MULT name = "CT"+NPL::itoa(i+1)+"_FRONT_RAW_MULT"; diff --git a/NPLib/Detectors/Samurai/TSamuraiFDC2Physics.cxx b/NPLib/Detectors/Samurai/TSamuraiFDC2Physics.cxx index 05818921b0505290870cee1e603c9fb9be1adead..74bbd25a62d25923031cad842af5f3aca1f96e08 100644 --- a/NPLib/Detectors/Samurai/TSamuraiFDC2Physics.cxx +++ b/NPLib/Detectors/Samurai/TSamuraiFDC2Physics.cxx @@ -1,5 +1,5 @@ /***************************************************************************** - * Copyright (C) 2009-2016 this file is part of the NPTool Project * + * Copyright (C) 2009-2020 this file is part of the NPTool Project * * * * For the licensing terms see $NPTOOL/Licence/NPTool_Licence * * For the list of contributors see $NPTOOL/Licence/Contributors * @@ -58,7 +58,7 @@ void TSamuraiFDC2Physics::BuildSimplePhysicalEvent(){ /////////////////////////////////////////////////////////////////////////// void TSamuraiFDC2Physics::BuildPhysicalEvent(){ PreTreat(); - //RemoveNoise(); + RemoveNoise(); // Map[detector & plane angle, vector of spatial information] static map<std::pair<unsigned int,double>, vector<double> > X ; @@ -68,7 +68,7 @@ void TSamuraiFDC2Physics::BuildPhysicalEvent(){ X.clear();Z.clear();R.clear(); unsigned int size = Detector.size(); for(unsigned int i = 0 ; i < size ; i++){ - if( DriftLength[i] > DriftLowThreshold && DriftLength[i] < DriftUpThreshold){ + if(DriftLength[i] > DriftLowThreshold && DriftLength[i] < DriftUpThreshold){ det = Detector[i]; layer = Layer[i]; wire = Wire[i]; @@ -88,9 +88,9 @@ void TSamuraiFDC2Physics::BuildPhysicalEvent(){ VX0.clear();VX100.clear(),D.clear(); for(auto it = X.begin();it!=X.end();++it){ D[it->first]=m_reconstruction.BuildTrack2D(X[it->first],Z[it->first],R[it->first],X0,X100,a,b); - std::ofstream f("distance.txt", std::ios::app); + /* std::ofstream f("distance.txt", std::ios::app); f<< D[it->first] << endl; - f.close(); + f.close();*/ // very large a means track perpendicular to the chamber, what happen when there is pile up if(abs(a)>1000) PileUp++; @@ -341,9 +341,9 @@ void TSamuraiFDC2Physics::AddDC(string name, NPL::XmlParser& xml){ else if(sDir=="Y") T=90*deg; else if(sDir=="U") - T=30*deg; - else if(sDir=="V") T=-30*deg; + else if(sDir=="V") + T=+30*deg; else{ cout << "ERROR: Unknown layer orientation for Samurai FDC2"<< endl; exit(1); diff --git a/NPLib/Detectors/Samurai/TSamuraiFDC2Physics.h b/NPLib/Detectors/Samurai/TSamuraiFDC2Physics.h index d7ffc98ba15c281fa28a0c21983ea91ab8f80dba..9259270c85daa85ef67afa9e5f612be2344cd440 100644 --- a/NPLib/Detectors/Samurai/TSamuraiFDC2Physics.h +++ b/NPLib/Detectors/Samurai/TSamuraiFDC2Physics.h @@ -1,7 +1,7 @@ #ifndef TSAMURAIFDC2PHYSICS_H #define TSAMURAIFDC2PHYSICS_H /***************************************************************************** - * Copyright (C) 2009-2016 this file is part of the NPTool Project * + * Copyright (C) 2009-2020 this file is part of the NPTool Project * * * * For the licensing terms see $NPTOOL/Licence/NPTool_Licence * * For the list of contributors see $NPTOOL/Licence/Contributors * diff --git a/NPLib/Detectors/eAGanil/TeAGanilPhysics.cxx b/NPLib/Detectors/eAGanil/TeAGanilPhysics.cxx index bb696e5bea270d8b3a427ee4bc85a178b2c4adde..bd55a402fda2c980d0bc4dfdd8eeeb9f283249e4 100644 --- a/NPLib/Detectors/eAGanil/TeAGanilPhysics.cxx +++ b/NPLib/Detectors/eAGanil/TeAGanilPhysics.cxx @@ -1,18 +1,18 @@ /***************************************************************************** - * Copyright (C) 2009-2020 this file is part of the NPTool Project * + * Copyright (C) 2009-2020 this file is part of the NPTool Project * * * * For the licensing terms see $NPTOOL/Licence/NPTool_Licence * * For the list of contributors see $NPTOOL/Licence/Contributors * *****************************************************************************/ /***************************************************************************** - * Original Author: Adrien Matta contact address: matta@lpccaen.in2p3.fr * + * Original Author: Adrien Matta contact address: matta@lpccaen.in2p3.fr * * * - * Creation Date : October 2020 * + * Creation Date : October 2020 * * Last update : * *---------------------------------------------------------------------------* * Decription: * - * This class hold eAGanil Treated data * + * This class hold eAGanil Treated data * * * *---------------------------------------------------------------------------* * Comment: * @@ -203,30 +203,18 @@ void TeAGanilPhysics::Clear() { /////////////////////////////////////////////////////////////////////////// void TeAGanilPhysics::ReadConfiguration(NPL::InputParser parser) { - vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("eAGanil"); + vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithTokenAndValue("eAGanil","Spectrometer"); if(NPOptionManager::getInstance()->GetVerboseLevel()) cout << "//// " << blocks.size() << " detectors found " << endl; - vector<string> cart = {"POS","Shape"}; - vector<string> sphe = {"R","Theta","Phi","Shape"}; - + vector<string> sphe = {"R","Theta","Phi","EntranceWidth","EntranceHeight","MomentumResolution"}; for(unsigned int i = 0 ; i < blocks.size() ; i++){ - if(blocks[i]->HasTokenList(cart)){ - if(NPOptionManager::getInstance()->GetVerboseLevel()) - cout << endl << "//// eAGanil " << i+1 << endl; - - TVector3 Pos = blocks[i]->GetTVector3("POS","mm"); - string Shape = blocks[i]->GetString("Shape"); - AddDetector(Pos,Shape); - } - else if(blocks[i]->HasTokenList(sphe)){ + if(blocks[i]->HasTokenList(sphe)){ if(NPOptionManager::getInstance()->GetVerboseLevel()) cout << endl << "//// eAGanil " << i+1 << endl; double R = blocks[i]->GetDouble("R","mm"); double Theta = blocks[i]->GetDouble("Theta","deg"); double Phi = blocks[i]->GetDouble("Phi","deg"); - string Shape = blocks[i]->GetString("Shape"); - AddDetector(R,Theta,Phi,Shape); } else{ cout << "ERROR: check your input file formatting " << endl; diff --git a/NPLib/TrackReconstruction/DCTest.cxx b/NPLib/TrackReconstruction/DCTest.cxx index c449ca20fb0ae7bc30cbe91e025aebdff0c0d9a5..dfb7a2ad63244f7cc68844a046a42bfd76f7d432 100644 --- a/NPLib/TrackReconstruction/DCTest.cxx +++ b/NPLib/TrackReconstruction/DCTest.cxx @@ -44,7 +44,7 @@ R__LOAD_LIBRARY(libNPTrackReconstruction.dylib) x.push_back(X); DrawCircle(X,Z,R); } - // add second track for pile up +/* // add second track for pile up a = rand.Uniform(-10,10); b = rand.Uniform(-100,100); DrawTrack(a,b,3); @@ -59,52 +59,41 @@ R__LOAD_LIBRARY(libNPTrackReconstruction.dylib) x.push_back(X); DrawCircle(X,Z,R); } - +*/ + // add noise + // Generate n points from the line +/* n = (int) rand.Uniform(1,3); + for(unsigned int i = 0 ; i < n ; i++){ + double Z = i*stepZ-100; + z.push_back(Z); + double R = rand.Uniform(-maxR,maxR); + r.push_back(R); + double X = rand.Uniform(-100,100); + x.push_back(X); + DrawCircle(X,Z,R); + } +*/ + double x0,x100,af,bf; double res = dcr.BuildTrack2D(x,z,r,x0,x100,af,bf); cout << res << endl; DrawTrack(af,bf,1); - // Check 2D track crossing point + // Check 2D Minimize Func c->cd(2); - bck->Draw(); - double ThetaL=-60*deg; - double ThetaH= 60*deg; - TVector3 L(rand.Uniform(-100,100),0,0); - TVector3 H(rand.Uniform(-100,100),0,0); - L.RotateZ(ThetaL); - H.RotateZ(ThetaH); - // direction of U and V wire - TVector3 u(0,1,0); - u.RotateZ(ThetaL); - - TVector3 v(0,1,0); - v.RotateZ(ThetaH); - - // Compute the coeff of the two line of vecotr u (v) going through H (L) - // dv : y = av*x+bv - double av = v.Y()/v.X(); - double bv = H.Y() - av*H.X(); - - DrawTrack(av,bv,1); - // du : y = au*x+bu - double au = u.Y()/u.X(); - double bu = L.Y() - au*L.X(); - - DrawTrack(au,bu,1); - TVector3 P; - dcr.ResolvePlane(L,ThetaL,H,ThetaH,P); - auto Lmk = new TMarker(L.X(),L.Y(),23); - auto Hmk = new TMarker(H.X(),H.Y(),27); - auto Pmk = new TMarker(P.X(),P.Y(),29); - Lmk->Draw(); - Hmk->Draw(); - Pmk->Draw(); + TGraph* ga = dcr.Scan(af,bf,0,af*0.5,af*2); + TGraph* gb = dcr.Scan(af,bf,1,bf*0.5,bf*2); + ga->SetLineColor(kAzure+7); + ga->Draw("ac"); + ga->GetYaxis()->SetRangeUser(0,1000); + gb->Draw("c"); + c->Update(); gPad->WaitPrimitive(); } } + // draw a track based on z=ax+b void DrawTrack(double a, double b, int style){ double x[2]={-100,100}; diff --git a/NPLib/TrackReconstruction/NPDCReconstruction.cxx b/NPLib/TrackReconstruction/NPDCReconstruction.cxx index 54e2efb463b16ad42df87c665ce1d140bcf81bc9..c1978c899dc4ecb6d9d5075cb5d762b831322c03 100644 --- a/NPLib/TrackReconstruction/NPDCReconstruction.cxx +++ b/NPLib/TrackReconstruction/NPDCReconstruction.cxx @@ -1,6 +1,3 @@ -#include"NPDCReconstruction.h" -#include "Math/Factory.h" -#include "TError.h" /***************************************************************************** * Copyright (C) 2009-2020 this file is part of the NPTool Project * * * @@ -19,6 +16,10 @@ * This class have all the method needed to analyse Drift Chambers * *****************************************************************************/ +#include"NPDCReconstruction.h" +#include "Math/Factory.h" +#include "TError.h" +#include "TGraph.h" using namespace std; using namespace NPL; @@ -30,16 +31,16 @@ DCReconstruction::DCReconstruction(){ m_min->SetPrintLevel(0); // this avoid error - gErrorIgnoreLevel = kError; + gErrorIgnoreLevel = kError; //m_min->SetMaxFunctionCalls(1000); //m_min->SetMaxIterations(1000); - //m_min->SetTolerance(1000); + //m_min->SetTolerance(1); //m_min->SetPrecision(1e-10); - } +} //////////////////////////////////////////////////////////////////////////////// DCReconstruction::~DCReconstruction(){ delete m_min; - } +} //////////////////////////////////////////////////////////////////////////////// double DCReconstruction::BuildTrack2D(const vector<double>& X,const vector<double>& Z,const vector<double>& R,double& X0,double& X100,double& a, double& b ){ @@ -47,23 +48,22 @@ double DCReconstruction::BuildTrack2D(const vector<double>& X,const vector<doubl fitZ=&Z; fitR=&R; // assume all X,Z,R of same size - size = X.size(); + sizeX = X.size(); // Define the starting point of the fit: a straight line passing through the // the first and last wire // z = ax+b -> x=(z-b)/a // ai = (Z[size-1]-Z[0])/(X[size-1]-R[size-1]-X[0]-R[0]); // bi = Z[0]-ai*(X[0]+R[0]); - ai = (Z[size-1]-Z[0])/(X[size-1]-X[0]); + ai = (Z[sizeX-1]-Z[0])/(X[sizeX-1]-X[0]); bi = Z[0]-ai*(X[0]); -// parameter[0]=ai; -// parameter[1]=bi; - m_min->Clear(); - m_min->SetLimitedVariable(0,"a",ai,1,ai*0.1,ai*10); - m_min->SetLimitedVariable(1,"b",bi,1,bi*0.1,bi*10); + m_min->SetVariable(0,"a",ai,1); + m_min->SetVariable(1,"b",bi,1); + // Perform minimisation m_min->Minimize(); + //std::cout << "EDM:" << m_min->Edm() << std::endl; // access set of parameter that produce the minimum const double *xs = m_min->X(); a=xs[0]; @@ -109,34 +109,54 @@ void DCReconstruction::ResolvePlane(const TVector3& L,const double& ThetaU ,cons xM=-10000; yM=-10000; } - PosXY.SetXYZ(xM,yM,0); - } //////////////////////////////////////////////////////////////////////////////// double DCReconstruction::SumD(const double* parameter ){ - unsigned int sizeX = fitX->size(); // Compute the sum P of the distance between the circle and the track P = 0; a = parameter[0]; b = parameter[1]; ab= a*b; a2=a*a; - for(unsigned int i = 0 ; i < sizeX ; i++){ c = (*fitX)[i]; d = (*fitZ)[i]; r = (*fitR)[i]; x = (a*d-ab+c)/(1+a2); z = a*x+b; - //P+= sqrt(abs( (x-c)*(x-c)+(z-d)*(z-d)-r*r)); - P+= abs( (x-c)*(x-c)+(z-d)*(z-d)-r*r); + p= (x-c)*(x-c)+(z-d)*(z-d)-r*r; + // numerical trick to have a smooth derivative instead of using abs + P+= sqrt(p*p+0.1); } - + // return normalized power - return P/size; + return P/sizeX; } +//////////////////////////////////////////////////////////////////////////////// +TGraph* DCReconstruction::Scan(double a, double b, int tovary, double minV, double maxV){ + vector<double> x,y; + unsigned int sizeT=1000; + double step = (maxV-minV)/sizeT; + double p[2]={a,b}; + for(unsigned int i = 0 ; i < sizeT ; i++){ + if(!tovary){ + p[0]=minV+step*i; + x.push_back(p[0]); + y.push_back(SumD(p)); + } + + else{ + p[1]=minV+step*i; + x.push_back(p[1]); + y.push_back(SumD(p)); + } + } + + TGraph* g = new TGraph(x.size(),&x[0],&y[0]); + return g; +} diff --git a/NPLib/TrackReconstruction/NPDCReconstruction.h b/NPLib/TrackReconstruction/NPDCReconstruction.h index 5659aaf2cbd3c30329472d3f779c45ee7e37f76d..4c478e9937ac94654ef107595c987c604b85562b 100644 --- a/NPLib/TrackReconstruction/NPDCReconstruction.h +++ b/NPLib/TrackReconstruction/NPDCReconstruction.h @@ -34,10 +34,11 @@ * intersaction in any given plane. This is done using ResolvePlane. * * - Resolving plane for two Z plane will provide a point and a direction * *****************************************************************************/ -#include<vector> -#include"TVector3.h" +#include <vector> +#include "TVector3.h" #include "Math/Minimizer.h" #include "Math/Functor.h" +class TGraph; namespace NPL{ class DCReconstruction{ @@ -57,6 +58,13 @@ namespace NPL{ // Function used by the minimizer in BuildTrack2D double SumD(const double* parameter ); + + // For debugging/optimisation + // Scan Sumd versus parameter a or b (tovary =0 for a, 1 for b) + // return a TGraph for display + TGraph* Scan(double a, double b, int tovary, double minV, double maxV); + + private: // private member used by SumD ROOT::Math::Minimizer* m_min; ROOT::Math::Functor m_func; @@ -64,8 +72,8 @@ namespace NPL{ const std::vector<double>* fitZ; const std::vector<double>* fitR; // used by SumD - unsigned int size ; - double P,a,b,ab,a2,c,d,x,z,r; + unsigned int sizeX ; + double P,p,a,b,ab,a2,c,d,x,z,r; // used by BuildTrack double ai,bi; double parameter[2]; @@ -73,7 +81,6 @@ namespace NPL{ long double av,bv,au,bu; double xM,yM; - }; } diff --git a/Projects/ComptonTelescope/online/src/online.cpp b/Projects/ComptonTelescope/online/src/online.cpp index 29a1dac79ac7dd81ba7da9d080e797cabe32d607..973c69703b9ae482901ab24d8761c01f7e88a049 100644 --- a/Projects/ComptonTelescope/online/src/online.cpp +++ b/Projects/ComptonTelescope/online/src/online.cpp @@ -17,7 +17,74 @@ #include <string> using namespace std; +//--//--// One-line setter for calorimeter //--//--// +void setCTCalorimeter(TComptonTelescopeData* ccamData, DecodeR* D, const int pixelNumber) +{ + ccamData -> SetCTCalorimeterTTowerNbr( 1 ); + ccamData -> SetCTCalorimeterTDetectorNbr( 1 );//Triggered ASIC number + ccamData -> SetCTCalorimeterTChannelNbr( D -> getPixelNumber() );//Pixel that triggered + ccamData -> SetCTCalorimeterTTime( D -> getTime() ); + for (int i = 0; i < pixelNumber; ++i) {//Loop on pixels + ccamData -> SetCTCalorimeterETowerNbr(1); + ccamData -> SetCTCalorimeterEDetectorNbr( 1 ); + ccamData -> SetCTCalorimeterEChannelNbr( i );//PMTÂ pixel number + ccamData -> SetCTCalorimeterEEnergy( D -> getData()[i] ); + }//End of loop on pixels +} + +//--//--// One-line setter for DSSSD(s) //--//--// + +// One-line setter for the Front of one DSSD +void setCTTrackerFront(TComptonTelescopeData* ccamData, newframe_t* event, int detNbr, int faceNbr, const int stripNumber) +{ + //detNbr and faceNbr are > 0 + ccamData -> SetCTTrackerFrontTTowerNbr(1); + ccamData -> SetCTTrackerFrontTDetectorNbr(detNbr); + ccamData -> SetCTTrackerFrontTStripNbr(33); + ccamData -> SetCTTrackerFrontTTime(event->timestamp); + for (int k = 0; k < stripNumber; k++) {//Loop on strips + ccamData -> SetCTTrackerFrontETowerNbr(1); + ccamData -> SetCTTrackerFrontEDetectorNbr(detNbr); + ccamData -> SetCTTrackerFrontEStripNbr(k+1); + ccamData -> SetCTTrackerFrontEEnergy(event->sample[detNbr-1][faceNbr-1][k]); + }//End of loop on strips +} + +// One-line setter for the Back of one DSSD +void setCTTrackerBack(TComptonTelescopeData* ccamData, newframe_t* event, int detNbr, int faceNbr, const int stripNumber) +{ + ccamData -> SetCTTrackerBackTTowerNbr(1); + ccamData -> SetCTTrackerBackTDetectorNbr(detNbr); + ccamData -> SetCTTrackerBackTStripNbr(33); + ccamData -> SetCTTrackerBackTTime(event->timestamp); + for (int k = 0; k < stripNumber; k++) {//Loop on strips + ccamData -> SetCTTrackerBackETowerNbr(1); + ccamData -> SetCTTrackerBackEDetectorNbr(detNbr); + ccamData -> SetCTTrackerBackEStripNbr(k+1); + ccamData -> SetCTTrackerBackEEnergy(event->sample[detNbr-1][faceNbr-1][k]); + }//End of loop on strips +} + +// One-line setter for DSSSD(s) +void setCTTracker(TComptonTelescopeData* ccamData, newframe_t* event, vector<int>* nb_asic, vector<int>* chain, const int stripNumber) +{ + for (vector<int>::iterator itchain = chain->begin(); itchain != chain->end(); ++itchain) {//Iterates on 2 faces + for (vector<int>::iterator itasic = nb_asic->begin(); itasic != nb_asic->end(); ++itasic) {//Iterates on 1 DSSSD + if (event->chip_data[*itchain][*itasic]) { // Test if data is present + switch (*itchain) { + case 0://Assuming 0 is front - to be checked + setCTTrackerFront(ccamData, event, *itasic+1, *itchain+1, stripNumber); + break; + case 1://Assuming 1 is back - to be checked + setCTTrackerBack(ccamData, event, *itasic+1, *itchain+1, stripNumber); + }//End switch + }//End if + }//End for + }//End for +} + +//--//--//--//--//--//--//--//--//--//--//--//--//--//--//--// int main() { /////////////////////////////////////////////////////////////////////////// @@ -49,48 +116,163 @@ int main() // read data file/flux and fill ccamData object std::cout << "Reading data\n"; - // instantiate DecodeR object reading calorimeter data flux - DecodeR* D = new DecodeR(false); + DecodeR* DR = new DecodeR(false); // Instantiates DecodeR object reading calorimeter data flux + DecodeD* DD = new DecodeD(true); // Instantiates DecodeD object reading DSSSD(s) data flux + newframe_t* event; - DecodeD* DD = new DecodeD(true); + //Sets where to look for data in DSSSD root frames + vector <int> chain ({0, 1});//Two faces + vector <int> nb_asic ({0});//First DSSSD - int i = 0; - int c = 0; + int i = 0;// ROSMAP files loop counter + int c = 0;// Event counter + // Set some constants const int pixelNumber = 64; - while (true) -// while (i<1) - { - // Load a file(s) - std::ifstream is; -// i = 0; - switch (i % 3) { - case 3: is.open("./mfm.bin", std::ios::binary); break; - case 4: is.open("./133Ba.bin", std::ios::binary); break; - case 5: is.open("./241Am.bin", std::ios::binary); break; - case 0: is.open("./207Bi.bin", std::ios::binary); break; - case 1: is.open("./241Am-1.bin", std::ios::binary); break; - case 2: is.open("./241Am-2.bin", std::ios::binary); break; - } - is.seekg (0, std::ios::end); - int length = is.tellg(); - is.seekg (0, ios::beg); - char* buffer = new char [length]; - is.read(buffer, length); - is.close(); - i++; - - // Read from file + const int stripNumber = 32; + const bool loopForever = true; + while (loopForever or i<1) + { + // Load a root file and setup DecodeD DD -> setTree("20200128_10h44_bi207_conv.root"); -// DD -> decodeEvent(); -// DD -> Dump(); int dlength = DD -> getLength(); - newframe_t* event; + + while (DD -> getCursor() < dlength and (loopForever or i<1)) + { + // Load a ROSMAP file + std::ifstream is; + i = 1; + switch (i % 3) { + case 3: is.open("./mfm.bin", std::ios::binary); break; + case 4: is.open("./133Ba.bin", std::ios::binary); break; + case 5: is.open("./241Am.bin", std::ios::binary); break; + case 0: is.open("./207Bi.bin", std::ios::binary); break; + case 1: is.open("./241Am-1.bin", std::ios::binary); break; + case 2: is.open("./241Am-2.bin", std::ios::binary); break; + } + is.seekg (0, std::ios::end); + int length = is.tellg(); + is.seekg (0, ios::beg); + char* buffer = new char [length]; + is.read(buffer, length); + is.close(); + i++; + + // Setup DecodeR + DR -> setRaw(buffer); + DR -> decodeRawMFM(); // get rid of the first two (empty) events + DR -> decodeRawMFM(); + + // Loop on some events + while (DD -> getCursor() < dlength and DR -> getCursor() < length) + { + // Clear raw and physics data + m_NPDetectorManager->ClearEventPhysics(); + m_NPDetectorManager->ClearEventData(); + + // Fill calorimeter data + DR -> decodeRawMFM(); + setCTCalorimeter(ccamData, DR, pixelNumber); + + // Fill DSSD data + DD -> decodeEvent(); + event = DD -> getEvent(); + setCTTracker(ccamData, event, &nb_asic, &chain, stripNumber); + + // Build physical event + m_NPDetectorManager->BuildPhysicalEvent(); + + // Fill object in output ROOT file + m_OutputTree->Fill(); + + // check spectra + m_NPDetectorManager->CheckSpectraServer(); + + c++; + //usleep(100);//Simulated 10kHz count rate + + }// End of loop on ROSMAP events + + delete [] buffer; + + }// End of loop on DSSSD data + + }// End of main loop + + // fill spectra + m_NPDetectorManager->WriteSpectra(); + + // delete Decoders + delete DR; + delete DD; + + // Essential + #if __cplusplus > 199711L && NPMULTITHREADING + m_NPDetectorManager->StopThread(); + #endif + RootOutput::Destroy(); + + return 0; +} + + + + + + +/* while (DD -> getCursor() < dlength) { + // Clear raw and physics data + m_NPDetectorManager->ClearEventPhysics(); + m_NPDetectorManager->ClearEventData(); + + //Read some data DD -> decodeEvent(); event = DD -> getEvent(); //Fill TComptonTelescopeData here (if possible) - } + for (vector<int>::iterator itchain = chain.begin(); itchain != chain.end(); ++itchain) {//Iterates on 2 faces + for (vector<int>::iterator itasic = nb_asic.begin(); itasic != nb_asic.end(); ++itasic) {//Iterates on 1 DSSSD + if (event->chip_data[*itchain][*itasic]) { // Test if data is present + switch (*itchain) { + case 0://Assuming 0 is front - to be checked + ccamData -> SetCTTrackerFrontTTowerNbr(1); + ccamData -> SetCTTrackerFrontTDetectorNbr(*itasic+1); + ccamData -> SetCTTrackerFrontTStripNbr(33); + ccamData -> SetCTTrackerFrontTTime(event->timestamp); + for (int k = 0; k < stripNumber; k++) {//Loop on strips + ccamData -> SetCTTrackerFrontETowerNbr(1); + ccamData -> SetCTTrackerFrontEDetectorNbr(*itasic+1); + ccamData -> SetCTTrackerFrontEStripNbr(k+1); + ccamData -> SetCTTrackerFrontEEnergy(event->sample[*itchain][*itasic][k]); + }//End of loop on strips + break; + case 1://Assuming 1 is back - to be checked + ccamData -> SetCTTrackerBackTTowerNbr(1); + ccamData -> SetCTTrackerBackTDetectorNbr(*itasic+1); + ccamData -> SetCTTrackerBackTStripNbr(33); + ccamData -> SetCTTrackerBackTTime(event->timestamp); + for (int k = 0; k < stripNumber; k++) {//Loop on strips + ccamData -> SetCTTrackerBackETowerNbr(1); + ccamData -> SetCTTrackerBackEDetectorNbr(*itasic+1); + ccamData -> SetCTTrackerBackEStripNbr(k+1); + ccamData -> SetCTTrackerBackEEnergy(event->sample[*itchain][*itasic][k]); + }//End of loop on strips + }//End switch + }//End if + }//End for + }//End for + + // Build physical event + m_NPDetectorManager->BuildPhysicalEvent(); + //Fill output tree + m_OutputTree->Fill(); + //check spectra + m_NPDetectorManager->CheckSpectraServer(); + + c++; + + }//End while + cout << "Read " << c << " DSSSD events from file" << endl; // Read from file D -> setRaw(buffer); @@ -132,26 +314,4 @@ int main() c++; usleep(100);//Simulated 10kHz count rate - } - - // delete buffer - delete [] buffer; - - } // end of main loop - - // fill spectra - m_NPDetectorManager->WriteSpectra(); - - // delete DecodeR - delete D; - - // Essential - #if __cplusplus > 199711L && NPMULTITHREADING - m_NPDetectorManager->StopThread(); - #endif - RootOutput::Destroy(); - - return 0; -} - - + }*/