diff --git a/NPLib/Detectors/Sofia/CMakeLists.txt b/NPLib/Detectors/Sofia/CMakeLists.txt index cf621e66f7216c94cc9e8ea5ec2962575e907bcd..65c96df0144859d82f0ae17a03d7f688b108d71c 100644 --- a/NPLib/Detectors/Sofia/CMakeLists.txt +++ b/NPLib/Detectors/Sofia/CMakeLists.txt @@ -5,6 +5,7 @@ add_custom_command(OUTPUT TSofTofWDataDict.cxx COMMAND ${CMAKE_BINARY_DIR}/scrip add_custom_command(OUTPUT TSofTofWPhysicsDict.cxx COMMAND ${CMAKE_BINARY_DIR}/scripts/build_dict.sh TSofTofWPhysics.h TSofTofWPhysicsDict.cxx TSofTofWPhysics.rootmap libNPSofia.dylib DEPENDS TSofTofWPhysics.h) add_custom_command(OUTPUT TSofMwpcDataDict.cxx COMMAND ${CMAKE_BINARY_DIR}/scripts/build_dict.sh TSofMwpcData.h TSofMwpcDataDict.cxx TSofMwpcData.rootmap libNPSofia.dylib DEPENDS TSofMwpcData.h) +add_custom_command(OUTPUT TSofMwpcPhysicsDict.cxx COMMAND ${CMAKE_BINARY_DIR}/scripts/build_dict.sh TSofMwpcPhysics.h TSofMwpcPhysicsDict.cxx TSofMwpcPhysics.rootmap libNPSofia.dylib DEPENDS TSofMwpcPhysics.h) add_custom_command(OUTPUT TSofAtDataDict.cxx COMMAND ${CMAKE_BINARY_DIR}/scripts/build_dict.sh TSofAtData.h TSofAtDataDict.cxx TSofAtData.rootmap libNPSofia.dylib DEPENDS TSofAtData.h) add_custom_command(OUTPUT TSofAtPhysicsDict.cxx COMMAND ${CMAKE_BINARY_DIR}/scripts/build_dict.sh TSofAtPhysics.h TSofAtPhysicsDict.cxx TSofAtPhysics.rootmap libNPSofia.dylib DEPENDS TSofAtPhysics.h) @@ -19,9 +20,9 @@ add_custom_command(OUTPUT TSofBeamIDDict.cxx COMMAND ${CMAKE_BINARY_DIR}/scripts add_custom_command(OUTPUT TSofFissionFragmentDict.cxx COMMAND ${CMAKE_BINARY_DIR}/scripts/build_dict.sh TSofFissionFragment.h TSofFissionFragmentDict.cxx TSofFissionFragment.rootmap libNPSofia.dylib DEPENDS TSofFissionFragment.h) -add_library(NPSofia SHARED TSofSciData.cxx TSofSciDataDict.cxx TSofSciPhysics.cxx TSofSciPhysicsDict.cxx TSofMwpcData.cxx TSofMwpcDataDict.cxx TSofAtData.cxx TSofAtDataDict.cxx TSofAtPhysics.cxx TSofAtPhysicsDict.cxx TSofTrimData.cxx TSofTrimDataDict.cxx TSofTrimPhysics.cxx TSofTrimPhysicsDict.cxx TSofTwimData.cxx TSofTwimDataDict.cxx TSofTwimPhysics.cxx TSofTwimPhysicsDict.cxx TSofTofWData.cxx TSofTofWDataDict.cxx TSofTofWPhysics.cxx TSofTofWPhysicsDict.cxx TSofBeamID.cxx TSofBeamIDDict.cxx TSofFissionFragment.cxx TSofFissionFragmentDict.cxx) +add_library(NPSofia SHARED TSofSciData.cxx TSofSciDataDict.cxx TSofSciPhysics.cxx TSofSciPhysicsDict.cxx TSofMwpcData.cxx TSofMwpcDataDict.cxx TSofMwpcPhysics.cxx TSofMwpcPhysicsDict.cxx TSofAtData.cxx TSofAtDataDict.cxx TSofAtPhysics.cxx TSofAtPhysicsDict.cxx TSofTrimData.cxx TSofTrimDataDict.cxx TSofTrimPhysics.cxx TSofTrimPhysicsDict.cxx TSofTwimData.cxx TSofTwimDataDict.cxx TSofTwimPhysics.cxx TSofTwimPhysicsDict.cxx TSofTofWData.cxx TSofTofWDataDict.cxx TSofTofWPhysics.cxx TSofTofWPhysicsDict.cxx TSofBeamID.cxx TSofBeamIDDict.cxx TSofFissionFragment.cxx TSofFissionFragmentDict.cxx) target_link_libraries(NPSofia ${ROOT_LIBRARIES} NPCore NPPhysics) -install(FILES TSofSciData.h TSofSciPhysics.h TSofMwpcData.h TSofAtData.h TSofAtPhysics.h TSofTrimData.h TSofTrimPhysics.h TSofTwimData.h TSofTwimPhysics.h TSofTofWData.h TSofTofWPhysics.h TSofBeamID.h TSofFissionFragment.h DESTINATION ${CMAKE_INCLUDE_OUTPUT_DIRECTORY}) +install(FILES TSofSciData.h TSofSciPhysics.h TSofMwpcData.h TSofMwpcPhysics.h TSofAtData.h TSofAtPhysics.h TSofTrimData.h TSofTrimPhysics.h TSofTwimData.h TSofTwimPhysics.h TSofTofWData.h TSofTofWPhysics.h TSofBeamID.h TSofFissionFragment.h DESTINATION ${CMAKE_INCLUDE_OUTPUT_DIRECTORY}) diff --git a/NPLib/Detectors/Sofia/TSofAtPhysics.cxx b/NPLib/Detectors/Sofia/TSofAtPhysics.cxx index db82c5a7872e401fb7ca0388c39ea8e902779ea1..749e6744833fa8ea02f35e8f4b1d04628e127634 100644 --- a/NPLib/Detectors/Sofia/TSofAtPhysics.cxx +++ b/NPLib/Detectors/Sofia/TSofAtPhysics.cxx @@ -80,19 +80,10 @@ void TSofAtPhysics::BuildPhysicalEvent() { PreTreat(); unsigned int mysizeE = m_PreTreatedData->GetMultiplicity(); - if(mysizeE != 4) - return; - double E[4]={-1,-1,-1,-1}; for (UShort_t e = 0; e < mysizeE ; e++) { - E[m_PreTreatedData->GetAnodeNbr(e)-1] = m_PreTreatedData->GetEnergy(e); - } - - if(E[0]>0 && E[1]>0 && E[2]>0 && E[3]>0){ - for(int i=0; i<4; i++){ - AnodeNbr.push_back(i+1); - Energy.push_back(E[i]); - } + AnodeNbr.push_back(m_PreTreatedData->GetAnodeNbr(e)); + Energy.push_back(m_PreTreatedData->GetEnergy(e)); } } diff --git a/NPLib/Detectors/Sofia/TSofBeamID.cxx b/NPLib/Detectors/Sofia/TSofBeamID.cxx index f078b93b0458c121e1df9fca58fc28f7cc0e4e56..2129233354dc4c3f8b27960e21dc048136446d9e 100644 --- a/NPLib/Detectors/Sofia/TSofBeamID.cxx +++ b/NPLib/Detectors/Sofia/TSofBeamID.cxx @@ -52,6 +52,7 @@ void TSofBeamID::Clear() { fBeam_Brho = -1; fBeam_XS2 = -1000; fBeam_XCC = -1000; + fBeam_YCC = -1000; } @@ -70,5 +71,6 @@ void TSofBeamID::Dump() const { cout << "fBeam_Brho: " << fBeam_Brho << endl; cout << "fBeam_XS2: " << fBeam_XS2 << endl; cout << "fBeam_XCC: " << fBeam_XCC << endl; + cout << "fBeam_YCC: " << fBeam_YCC << endl; } diff --git a/NPLib/Detectors/Sofia/TSofBeamID.h b/NPLib/Detectors/Sofia/TSofBeamID.h index be129d80698ecc9cdcfe0fd7c958943fd0c5e532..ba810d5331781712b1ca4b787d61e33605792bf2 100644 --- a/NPLib/Detectors/Sofia/TSofBeamID.h +++ b/NPLib/Detectors/Sofia/TSofBeamID.h @@ -43,6 +43,7 @@ class TSofBeamID : public TObject { double fBeam_Brho; double fBeam_XS2; double fBeam_XCC; + double fBeam_YCC; ////////////////////////////////////////////////////////////// // Constructor and destructor @@ -75,6 +76,7 @@ class TSofBeamID : public TObject { inline void SetBrho(double val){fBeam_Brho = val;};//! inline void SetXS2(double val){fBeam_XS2 = val;};//! inline void SetXCC(double val){fBeam_XCC = val;};//! + inline void SetYCC(double val){fBeam_YCC = val;};//! ////////////////////// GETTERS //////////////////////// inline double GetZbeam() const {return fBeam_Z;}//! @@ -86,6 +88,7 @@ class TSofBeamID : public TObject { inline double GetBrho() const {return fBeam_Brho;}//! inline double GetXS2() const {return fBeam_XS2;}//! inline double GetXCC() const {return fBeam_XCC;}//! + inline double GetYCC() const {return fBeam_YCC;}//! ////////////////////////////////////////////////////////////// // Required for ROOT dictionnary diff --git a/NPLib/Detectors/Sofia/TSofMwpcPhysics.cxx b/NPLib/Detectors/Sofia/TSofMwpcPhysics.cxx new file mode 100644 index 0000000000000000000000000000000000000000..182332b55ed56692c1044605c8d6b3944bfbf8b5 --- /dev/null +++ b/NPLib/Detectors/Sofia/TSofMwpcPhysics.cxx @@ -0,0 +1,380 @@ +/***************************************************************************** + * 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: Pierre Morfouace contact address: pierre.morfouace2@cea.fr * + * * + * Creation Date : November 2020 * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * This class hold SofMwpc Treated data * + * * + *---------------------------------------------------------------------------* + * Comment: * + * * + * * + *****************************************************************************/ + +#include "TSofMwpcPhysics.h" + +// STL +#include <sstream> +#include <iostream> +#include <cmath> +#include <stdlib.h> +#include <limits> +using namespace std; + +// NPL +#include "RootInput.h" +#include "RootOutput.h" +#include "NPDetectorFactory.h" +#include "NPOptionManager.h" + +// ROOT +#include "TChain.h" + +ClassImp(TSofMwpcPhysics) + + + /////////////////////////////////////////////////////////////////////////// +TSofMwpcPhysics::TSofMwpcPhysics() + : m_EventData(new TSofMwpcData), + m_PreTreatedData(new TSofMwpcData), + m_EventPhysics(this), + m_NumberOfDetectors(0){ + } + +/////////////////////////////////////////////////////////////////////////// +/// A usefull method to bundle all operation to add a detector +void TSofMwpcPhysics::AddDetector(TVector3 Pos){ + // In That simple case nothing is done + // Typically for more complex detector one would calculate the relevant + // positions (stripped silicon) or angles (gamma array) + DetPosX.push_back(Pos.X()); + DetPosY.push_back(Pos.Y()); + DetPosZ.push_back(Pos.Z()); + + m_NumberOfDetectors++; +} + +/////////////////////////////////////////////////////////////////////////// +void TSofMwpcPhysics::AddDetector(double R, double Theta, double Phi){ + // Compute the TVector3 corresponding + TVector3 Pos(R*sin(Theta)*cos(Phi),R*sin(Theta)*sin(Phi),R*cos(Theta)); + // Call the cartesian method + AddDetector(Pos); +} + +/////////////////////////////////////////////////////////////////////////// +void TSofMwpcPhysics::BuildSimplePhysicalEvent() { + BuildPhysicalEvent(); +} + + + +/////////////////////////////////////////////////////////////////////////// +void TSofMwpcPhysics::BuildPhysicalEvent() { + + PreTreat(); + + unsigned int mysizeE = m_PreTreatedData->GetMultiplicity(); + + vector<double> ChargeArray; + ChargeArray.resize(288,0); + + int StripMaxX1[4]={0,0,0,0}; + int StripMaxX2[4]={0,0,0,0}; + int StripMaxY[4] ={0,0,0,0}; + double ChargeMaxX1[4]={0,0,0,0}; + double ChargeMaxX2[4]={0,0,0,0}; + double ChargeMaxY[4] ={0,0,0,0}; + + + for(int i=0; i<m_NumberOfDetectors; i++){ + Buffer_X1_Q.push_back(ChargeArray); + Buffer_X2_Q.push_back(ChargeArray); + Buffer_Y_Q.push_back(ChargeArray); + } + + for (UShort_t e = 0; e < mysizeE ; e++) { + int det = m_PreTreatedData->GetDetectorNbr(e); + int plane = m_PreTreatedData->GetPlane(e); + int strip = m_PreTreatedData->GetPad(e);// starts at 0 + double charge = m_PreTreatedData->GetCharge(e); + + // Xup for MWPC2 and 3 and X for MWPC1 and 4 + if(plane==1){ + Buffer_X1_Q[det-1][strip] = charge; + + if(charge>ChargeMaxX1[det-1]){ + ChargeMaxX1[det-1] = charge; + StripMaxX1[det-1] = strip; + } + } + // Xdown for MWPC2 and 3 and nothing for MWPC1 and 4 + else if(plane==2){ + Buffer_X2_Q[det-1][strip] = charge; + + if(charge>ChargeMaxX2[det-1]){ + ChargeMaxX2[det-1] = charge; + StripMaxX2[det-1] = strip; + } + } + // Y for all MWPCx + if(plane==3){ + Buffer_Y_Q[det-1][strip] = charge; + + if(charge>ChargeMaxY[det-1]){ + ChargeMaxY[det-1] = charge; + StripMaxY[det-1] = strip; + } + } + } + + double X1 = -100; + double X2 = -100; + double Y = -100; + for(int i=0; i<m_NumberOfDetectors; i++){ + double qleft = Buffer_X1_Q[i][StripMaxX1[i]-1]; + double qright = Buffer_X1_Q[i][StripMaxX1[i]+1]; + if(qleft>0 && qright>0){ + X1 = GetPositionX(ChargeMaxX1[i], StripMaxX1[i], qleft, qright); + X1 = X1 - DetPosX[i]; + } + qleft = Buffer_X2_Q[i][StripMaxX2[i]-1]; + qright = Buffer_X2_Q[i][StripMaxX2[i]+1]; + if(qleft>0 && qright>0){ + X2 = GetPositionX(ChargeMaxX2[i], StripMaxX2[i], qleft, qright); + X2 = X2 - DetPosX[i]; + } + double qdown = Buffer_Y_Q[i][StripMaxY[i]-1]; + double qup = Buffer_Y_Q[i][StripMaxY[i]+1]; + if(qdown>0 && qup>0){ + Y = GetPositionY(ChargeMaxY[i], StripMaxY[i], qdown, qup); + Y = Y - DetPosY[i]; + } + + DetectorNbr.push_back(i+1); + PositionX1.push_back(X1); + PositionX2.push_back(X2); + PositionY.push_back(Y); + } +} + +/////////////////////////////////////////////////////////////////////////// +double TSofMwpcPhysics::GetPositionX(double qmax, int padmax, double qleft, double qright){ + double fwx = 3.125; + double fSize = 200.0; + + double a3 = TMath::Pi() * fwx / (TMath::ACosH(0.5 * (TMath::Sqrt(qmax / qleft) + TMath::Sqrt(qmax / qright)))); + + Double_t a2 = (a3 / TMath::Pi()) * TMath::ATanH((TMath::Sqrt(qmax / qleft) - TMath::Sqrt(qmax / qright)) / + (2 * TMath::SinH(TMath::Pi() * fwx / a3))); + + return (-1. * padmax * fwx + (fSize / 2) - (fwx / 2) - a2); // Left is positive and right negative +} + +/////////////////////////////////////////////////////////////////////////// +double TSofMwpcPhysics::GetPositionY(double qmax, int padmax, double qdown, double qup){ + double fwy = 3.125; + double fSize = 200.0; + double a3 = TMath::Pi() * fwy / (TMath::ACosH(0.5 * (TMath::Sqrt(qmax / qdown) + TMath::Sqrt(qmax / qup)))); + + Double_t a2 = (a3 / TMath::Pi()) * TMath::ATanH((TMath::Sqrt(qmax / qdown) - TMath::Sqrt(qmax / qup)) / + (2 * TMath::SinH(TMath::Pi() * fwy / a3))); + + return (padmax * fwy - (fSize / 2) + (fwy / 2) + a2); +} + +/////////////////////////////////////////////////////////////////////////// +void TSofMwpcPhysics::PreTreat() { + // This method typically applies thresholds and calibrations + // Might test for disabled channels for more complex detector + + // clear pre-treated object + ClearPreTreatedData(); + + // instantiate CalibrationManager + static CalibrationManager* Cal = CalibrationManager::getInstance(); + + unsigned int mysize = m_EventData->GetMultiplicity(); + for (unsigned int i = 0; i < mysize ; ++i) { + m_PreTreatedData->SetDetectorNbr(m_EventData->GetDetectorNbr(i)); + m_PreTreatedData->SetPlane(m_EventData->GetPlane(i)); + m_PreTreatedData->SetPad(m_EventData->GetPad(i)); + m_PreTreatedData->SetCharge(m_EventData->GetCharge(i)); + } +} + + + +/////////////////////////////////////////////////////////////////////////// +void TSofMwpcPhysics::ReadAnalysisConfig() { + bool ReadingStatus = false; + + // path to file + string FileName = "./configs/ConfigSofMwpc.dat"; + + // open analysis config file + ifstream AnalysisConfigFile; + AnalysisConfigFile.open(FileName.c_str()); + + if (!AnalysisConfigFile.is_open()) { + cout << " No ConfigSofMwpc.dat found: Default parameter loaded for Analayis " << FileName << endl; + return; + } + cout << " Loading user parameter for Analysis from ConfigSofMwpc.dat " << endl; + + // Save it in a TAsciiFile + TAsciiFile* asciiConfig = RootOutput::getInstance()->GetAsciiFileAnalysisConfig(); + asciiConfig->AppendLine("%%% ConfigSofMwpc.dat %%%"); + asciiConfig->Append(FileName.c_str()); + asciiConfig->AppendLine(""); + // read analysis config file + string LineBuffer,DataBuffer,whatToDo; + while (!AnalysisConfigFile.eof()) { + // Pick-up next line + getline(AnalysisConfigFile, LineBuffer); + + // search for "header" + string name = "ConfigSofMwpc"; + if (LineBuffer.compare(0, name.length(), name) == 0) + ReadingStatus = true; + + // loop on tokens and data + while (ReadingStatus ) { + whatToDo=""; + AnalysisConfigFile >> whatToDo; + + // Search for comment symbol (%) + if (whatToDo.compare(0, 1, "%") == 0) { + AnalysisConfigFile.ignore(numeric_limits<streamsize>::max(), '\n' ); + } + + else if (whatToDo=="E_THRESHOLD") { + AnalysisConfigFile >> DataBuffer; + m_E_Threshold = atof(DataBuffer.c_str()); + cout << whatToDo << " " << m_E_Threshold << endl; + } + + else { + ReadingStatus = false; + } + } + } +} + + +/////////////////////////////////////////////////////////////////////////// +void TSofMwpcPhysics::Clear() { + DetectorNbr.clear(); + PositionX1.clear(); + PositionX2.clear(); + PositionY.clear(); + + Buffer_X1_Q.clear(); + Buffer_X2_Q.clear(); + Buffer_Y_Q.clear(); +} + + + +/////////////////////////////////////////////////////////////////////////// +void TSofMwpcPhysics::ReadConfiguration(NPL::InputParser parser) { + vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("SofMwpc"); + if(NPOptionManager::getInstance()->GetVerboseLevel()) + cout << "//// " << blocks.size() << " detectors found " << endl; + + vector<string> cart = {"POS"}; + vector<string> sphe = {"R","Theta","Phi"}; + + for(unsigned int i = 0 ; i < blocks.size() ; i++){ + if(blocks[i]->HasTokenList(cart)){ + if(NPOptionManager::getInstance()->GetVerboseLevel()) + cout << endl << "//// SofMwpc " << i+1 << endl; + + TVector3 Pos = blocks[i]->GetTVector3("POS","mm"); + AddDetector(Pos); + } + else if(blocks[i]->HasTokenList(sphe)){ + if(NPOptionManager::getInstance()->GetVerboseLevel()) + cout << endl << "//// SofMwpc " << i+1 << endl; + double R = blocks[i]->GetDouble("R","mm"); + double Theta = blocks[i]->GetDouble("Theta","deg"); + double Phi = blocks[i]->GetDouble("Phi","deg"); + AddDetector(R,Theta,Phi); + } + else{ + cout << "ERROR: check your input file formatting " << endl; + exit(1); + } + } + + ReadAnalysisConfig(); +} + + +/////////////////////////////////////////////////////////////////////////// +void TSofMwpcPhysics::AddParameterToCalibrationManager() { + CalibrationManager* Cal = CalibrationManager::getInstance(); + + for(int sec = 0; sec < m_NumberOfDetectors; sec++){ + Cal->AddParameter("SofMwpc","SEC"+NPL::itoa(sec+1)+"_ALIGN","SofMwpc_SEC"+NPL::itoa(sec+1)+"_ALIGN"); + } +} + +/////////////////////////////////////////////////////////////////////////// +void TSofMwpcPhysics::InitializeRootInputRaw() { + TChain* inputChain = RootInput::getInstance()->GetChain(); + inputChain->SetBranchStatus("SofMwpc", true ); + inputChain->SetBranchAddress("SofMwpc", &m_EventData ); +} + + + +/////////////////////////////////////////////////////////////////////////// +void TSofMwpcPhysics::InitializeRootInputPhysics() { + TChain* inputChain = RootInput::getInstance()->GetChain(); + inputChain->SetBranchAddress("SofMwpc", &m_EventPhysics); +} + + + +/////////////////////////////////////////////////////////////////////////// +void TSofMwpcPhysics::InitializeRootOutput() { + TTree* outputTree = RootOutput::getInstance()->GetTree(); + outputTree->Branch("SofMwpc", "TSofMwpcPhysics", &m_EventPhysics); +} + + + +//////////////////////////////////////////////////////////////////////////////// +// Construct Method to be pass to the DetectorFactory // +//////////////////////////////////////////////////////////////////////////////// +NPL::VDetector* TSofMwpcPhysics::Construct() { + return (NPL::VDetector*) new TSofMwpcPhysics(); +} + + + +//////////////////////////////////////////////////////////////////////////////// +// Registering the construct method to the factory // +//////////////////////////////////////////////////////////////////////////////// +extern "C"{ + class proxy_SofMwpc{ + public: + proxy_SofMwpc(){ + NPL::DetectorFactory::getInstance()->AddToken("SofMwpc","SofMwpc"); + NPL::DetectorFactory::getInstance()->AddDetector("SofMwpc",TSofMwpcPhysics::Construct); + } + }; + + proxy_SofMwpc p_SofMwpc; +} + diff --git a/NPLib/Detectors/Sofia/TSofMwpcPhysics.h b/NPLib/Detectors/Sofia/TSofMwpcPhysics.h new file mode 100644 index 0000000000000000000000000000000000000000..14935bd46a5c3700a0792f68bda9b1e301308c73 --- /dev/null +++ b/NPLib/Detectors/Sofia/TSofMwpcPhysics.h @@ -0,0 +1,162 @@ +#ifndef TSofMwpcPHYSICS_H +#define TSofMwpcPHYSICS_H +/***************************************************************************** + * 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: Pierre Morfouace contact address: pierre.morfouace2@cea.fr * + * * + * Creation Date : November 2020 * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * This class hold SofMwpc Treated data * + * * + *---------------------------------------------------------------------------* + * Comment: * + * * + * * + *****************************************************************************/ + +// C++ headers +#include <vector> +#include <map> +#include <string> +using namespace std; + +// ROOT headers +#include "TObject.h" +#include "TH1.h" +#include "TVector3.h" +#include "TSpline.h" +// NPTool headers +#include "TSofMwpcData.h" +#include "NPCalibrationManager.h" +#include "NPVDetector.h" +#include "NPInputParser.h" + + + +class TSofMwpcPhysics : public TObject, public NPL::VDetector { + ////////////////////////////////////////////////////////////// + // constructor and destructor + public: + TSofMwpcPhysics(); + ~TSofMwpcPhysics() {}; + + + ////////////////////////////////////////////////////////////// + // Inherited from TObject and overriden to avoid warnings + public: + void Clear(); + void Clear(const Option_t*) {}; + + + ////////////////////////////////////////////////////////////// + // data obtained after BuildPhysicalEvent() and stored in + // output ROOT file + public: + vector<int> DetectorNbr; + vector<double> PositionX1; + vector<double> PositionX2; + vector<double> PositionY; + + private: + vector<double> DetPosX;//! + vector<double> DetPosY;//! + vector<double> DetPosZ;//! + vector<vector<double>> Buffer_X1_Q;//! + vector<vector<double>> Buffer_X2_Q;//! + vector<vector<double>> Buffer_Y_Q;//! + + /// A usefull method to bundle all operation to add a detector + void AddDetector(TVector3 Pos); + void AddDetector(double R, double Theta, double Phi); + + ////////////////////////////////////////////////////////////// + // methods inherited from the VDetector ABC class + public: + // read stream from ConfigFile to pick-up detector parameters + void ReadConfiguration(NPL::InputParser); + + // add parameters to the CalibrationManger + void AddParameterToCalibrationManager(); + + // method called event by event, aiming at extracting the + // physical information from detector + void BuildPhysicalEvent(); + + // same as BuildPhysicalEvent() method but with a simpler + // treatment + void BuildSimplePhysicalEvent(); + + // same as above but for online analysis + void BuildOnlinePhysicalEvent() {BuildPhysicalEvent();}; + + // activate raw data object and branches from input TChain + // in this method mother branches (Detector) AND daughter leaves + // (fDetector_parameter) have to be activated + void InitializeRootInputRaw(); + + // activate physics data object and branches from input TChain + // in this method mother branches (Detector) AND daughter leaves + // (fDetector_parameter) have to be activated + void InitializeRootInputPhysics(); + + // create branches of output ROOT file + void InitializeRootOutput(); + + // clear the raw and physical data objects event by event + void ClearEventPhysics() {Clear();} + void ClearEventData() {m_EventData->Clear();} + + + ////////////////////////////////////////////////////////////// + // specific methods to SofMwpc array + public: + // remove bad channels, calibrate the data and apply thresholds + void PreTreat(); + + // clear the pre-treated object + void ClearPreTreatedData() {m_PreTreatedData->Clear();} + + // read the user configuration file. If no file is found, load standard one + void ReadAnalysisConfig(); + + // give and external TSofMwpcData object to TSofMwpcPhysics. + // needed for online analysis for example + void SetRawDataPointer(TSofMwpcData* rawDataPointer) {m_EventData = rawDataPointer;} + + double GetPositionX(double qmax, int padmax, double qleft, double qright); + double GetPositionY(double qmax, int padmax, double qdown, double qup); + + // objects are not written in the TTree + private: + TSofMwpcData* m_EventData; //! + TSofMwpcData* m_PreTreatedData; //! + TSofMwpcPhysics* m_EventPhysics; //! + + // getters for raw and pre-treated data object + public: + TSofMwpcData* GetRawData() const {return m_EventData;} + TSofMwpcData* GetPreTreatedData() const {return m_PreTreatedData;} + + // parameters used in the analysis + private: + double m_E_Threshold; //! + + // number of detectors + private: + int m_NumberOfDetectors; //! + + // Static constructor to be passed to the Detector Factory + public: + static NPL::VDetector* Construct(); + + ClassDef(TSofMwpcPhysics,1) // SofMwpcPhysics structure +}; +#endif diff --git a/NPLib/Detectors/Sofia/TSofTrimPhysics.cxx b/NPLib/Detectors/Sofia/TSofTrimPhysics.cxx index 35a50d1722b53265f04e027b2d746c495571cba7..b2f606d3439d7a4a0059f1911f931b18d168117a 100644 --- a/NPLib/Detectors/Sofia/TSofTrimPhysics.cxx +++ b/NPLib/Detectors/Sofia/TSofTrimPhysics.cxx @@ -48,7 +48,9 @@ TSofTrimPhysics::TSofTrimPhysics() m_PreTreatedData(new TSofTrimData), m_EventPhysics(this), m_NumberOfDetectors(0), - m_Beta(-1), + m_Beta(-1), + m_IsSplineSectionDriftTime(false), + m_IsSplineSectionAngle(false), m_BetaNorm(0.838), m_NumberOfSections(3), m_NumberOfAnodesPaired(3), @@ -194,7 +196,7 @@ void TSofTrimPhysics::BuildPhysicalEvent() { Ep3[i] = norm3 * Ep3[i] / (p0_3[i] + p1_3[i]*TMath::Power(m_Beta, -5./3.)); // Angle correction per pair: spline - Ep1[i] = Ep1[i] / fcorr_EvsA[i][0]->Eval(Ddt) * fcorr_EvsA[i][0]->Eval(0); + /*Ep1[i] = Ep1[i] / fcorr_EvsA[i][0]->Eval(Ddt) * fcorr_EvsA[i][0]->Eval(0); Ep2[i] = Ep2[i] / fcorr_EvsA[i][1]->Eval(Ddt) * fcorr_EvsA[i][1]->Eval(0); Ep3[i] = Ep3[i] / fcorr_EvsA[i][2]->Eval(Ddt) * fcorr_EvsA[i][2]->Eval(0); @@ -202,12 +204,17 @@ void TSofTrimPhysics::BuildPhysicalEvent() { Ep1[i] = Ep1[i] / fcorr_EvsDT[i][0]->Eval(DTp1[i]) * fcorr_EvsDT[i][0]->Eval(3000); Ep2[i] = Ep2[i] / fcorr_EvsDT[i][1]->Eval(DTp2[i]) * fcorr_EvsDT[i][1]->Eval(3000); Ep3[i] = Ep3[i] / fcorr_EvsDT[i][2]->Eval(DTp3[i]) * fcorr_EvsDT[i][2]->Eval(3000); - +*/ // Summing up Anode Energy per section Esec[i] = (Ep1[i] + Ep2[i] + Ep3[i])/3; + // Angle correction per section: spline + if(m_IsSplineSectionAngle) + Esec[i] = Esec[i] / fcorr_sec_angle[i]->Eval(Ddt) * fcorr_sec_angle[i]->Eval(0); + // 2nd DT correction per section: spline - Esec[i] = Esec[i] / fcorr_sec[i]->Eval(DTp2[i]) * fcorr_sec[i]->Eval(0); + if(m_IsSplineSectionDriftTime) + Esec[i] = Esec[i] / fcorr_sec_dt[i]->Eval(DTp2[i]) * fcorr_sec_dt[i]->Eval(0); // Section ALignement Esec[i] = Cal->ApplyCalibration("SofTrim/SEC"+NPL::itoa(i+1)+"_ALIGN",Esec[i]); @@ -308,17 +315,25 @@ void TSofTrimPhysics::ReadAnalysisConfig() { else if (whatToDo=="SPLINE_PAIR_DT_PATH") { AnalysisConfigFile >> DataBuffer; m_SPLINE_PAIR_DT_PATH = DataBuffer; - cout << "*** Loading Spline for Dritf Time correction per pair ***" << endl; + cout << "*** Loading Spline for Drtft Time correction per pair ***" << endl; LoadSplinePairDriftTime(); } else if (whatToDo=="SPLINE_SECTION_DT_PATH") { AnalysisConfigFile >> DataBuffer; m_SPLINE_SECTION_DT_PATH = DataBuffer; - cout << "*** Loading Spline for Dritf Time correction per section ***" << endl; + cout << "*** Loading Spline for Drift Time correction per section ***" << endl; LoadSplineSectionDriftTime(); } + else if (whatToDo=="SPLINE_SECTION_ANGLE_PATH") { + AnalysisConfigFile >> DataBuffer; + m_SPLINE_SECTION_ANGLE_PATH = DataBuffer; + cout << "*** Loading Spline for Angle correction per section ***" << endl; + LoadSplineSectionAngle(); + } + + else if (whatToDo=="E_THRESHOLD") { AnalysisConfigFile >> DataBuffer; @@ -389,22 +404,42 @@ void TSofTrimPhysics::LoadSplinePairDriftTime(){ /////////////////////////////////////////////////////////////////////////// void TSofTrimPhysics::LoadSplineSectionDriftTime(){ - TString filename = m_SPLINE_SECTION_DT_PATH; - TFile* ifile = new TFile(filename,"read"); + TString filename = m_SPLINE_SECTION_DT_PATH; + TFile* ifile = new TFile(filename,"read"); if(ifile->IsOpen()){ cout << "Loading splines..." << endl; for(int s=0; s<m_NumberOfSections; s++){ TString splinename = Form("spline_sec%i",s+1); - fcorr_sec[s] = (TSpline3*) ifile->FindObjectAny(splinename); - cout << fcorr_sec[s]->GetName() << endl; + fcorr_sec_dt[s] = (TSpline3*) ifile->FindObjectAny(splinename); + cout << fcorr_sec_dt[s]->GetName() << endl; } + m_IsSplineSectionDriftTime = true; } else cout << "File " << filename << " not found!" << endl; ifile->Close(); } +/////////////////////////////////////////////////////////////////////////// +void TSofTrimPhysics::LoadSplineSectionAngle(){ + TString filename = m_SPLINE_SECTION_ANGLE_PATH; + TFile* ifile = new TFile(filename,"read"); + + if(ifile->IsOpen()){ + cout << "Loading splines..." << endl; + for(int s=0; s<m_NumberOfSections; s++){ + TString splinename = Form("spline_sec%i",s+1); + + fcorr_sec_angle[s] = (TSpline3*) ifile->FindObjectAny(splinename); + cout << fcorr_sec_angle[s]->GetName() << endl; + } + m_IsSplineSectionAngle = true; + } + else + cout << "File " << filename << " not found!" << endl; + ifile->Close(); +} /////////////////////////////////////////////////////////////////////////// void TSofTrimPhysics::Clear() { diff --git a/NPLib/Detectors/Sofia/TSofTrimPhysics.h b/NPLib/Detectors/Sofia/TSofTrimPhysics.h index 923b46288ada001d9abb38d94f3c266a76eeaade..7053b20fd4188c2ddbdd6d471f44a4dbe4bbb10e 100644 --- a/NPLib/Detectors/Sofia/TSofTrimPhysics.h +++ b/NPLib/Detectors/Sofia/TSofTrimPhysics.h @@ -131,6 +131,7 @@ class TSofTrimPhysics : public TObject, public NPL::VDetector { void LoadSplinePairAngle(); void LoadSplinePairDriftTime(); void LoadSplineSectionDriftTime(); + void LoadSplineSectionAngle(); void SetBeta(double beta) {m_Beta = beta;} double GetBeta() {return m_Beta;} @@ -156,11 +157,15 @@ class TSofTrimPhysics : public TObject, public NPL::VDetector { string m_SPLINE_PAIR_ANGLE_PATH; //! string m_SPLINE_PAIR_DT_PATH; //! string m_SPLINE_SECTION_DT_PATH; //! + string m_SPLINE_SECTION_ANGLE_PATH; //! + bool m_IsSplineSectionDriftTime; //! + bool m_IsSplineSectionAngle; //! private: TSpline3* fcorr_EvsA[3][3]; //! TSpline3* fcorr_EvsDT[3][3]; //! - TSpline3* fcorr_sec[3]; //! + TSpline3* fcorr_sec_angle[3]; //! + TSpline3* fcorr_sec_dt[3]; //! // number of detectors private: