From e5446d734fa5c05f1b216001f570c32ad1d7778f Mon Sep 17 00:00:00 2001
From: Adrien Matta <matta@lpccaen.in2p3.fr>
Date: Tue, 31 May 2022 11:14:22 +0200
Subject: [PATCH] * Fixing SuperX3 placement by Spherical coordinate
---
NPLib/Core/RootInput.cxx | 202 +++++++-------
NPLib/Detectors/MDM/CMakeLists.txt | 2 +-
NPLib/Detectors/Minos/CMakeLists.txt | 2 +-
NPLib/Detectors/Minos/TMinosPhysics.cxx | 248 +++++++++---------
NPLib/TrackReconstruction/CMakeLists.txt | 2 +-
.../TrackReconstruction/NPLinearRansac3D.cxx | 234 +++++++++--------
NPSimulation/Detectors/SuperX3/SuperX3.cc | 16 +-
Projects/SuperX3/SX3_20cm.detector | 26 ++
8 files changed, 377 insertions(+), 355 deletions(-)
create mode 100644 Projects/SuperX3/SX3_20cm.detector
diff --git a/NPLib/Core/RootInput.cxx b/NPLib/Core/RootInput.cxx
index b92a4f766..a10a49cce 100644
--- a/NPLib/Core/RootInput.cxx
+++ b/NPLib/Core/RootInput.cxx
@@ -19,22 +19,22 @@
* *
*****************************************************************************/
-#include <iostream>
-#include <sstream>
+#include <dlfcn.h>
#include <fstream>
+#include <iostream>
#include <limits>
-#include <sys/stat.h>
+#include <sstream>
#include <stdlib.h>
-#include <dlfcn.h>
+#include <sys/stat.h>
+#include "NPOptionManager.h"
#include "RootInput.h"
#include "TAsciiFile.h"
-#include "NPOptionManager.h"
using namespace std;
RootInput* RootInput::instance = 0;
////////////////////////////////////////////////////////////////////////////////
-RootInput* RootInput::getInstance(std::string configFileName){
+RootInput* RootInput::getInstance(std::string configFileName) {
// A new instance of RootInput is created if it does not exist:
if (instance == 0) {
instance = new RootInput(configFileName);
@@ -45,7 +45,7 @@ RootInput* RootInput::getInstance(std::string configFileName){
}
////////////////////////////////////////////////////////////////////////////////
-void RootInput::Destroy(){
+void RootInput::Destroy() {
if (instance != 0) {
delete instance;
instance = 0;
@@ -53,64 +53,64 @@ void RootInput::Destroy(){
}
////////////////////////////////////////////////////////////////////////////////
-void RootInput::ReadOldStyleInputFile(NPL::InputParser& parser){
- vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("TTreeName");
- pTreeName = blocks[0]->GetLines()[0];
- vector<NPL::InputBlock*> files = parser.GetAllBlocksWithToken("RootFileName");
- if(files.size()>0){
- vector<string> lines=files[0]->GetLines();
- unsigned int size = lines.size();
- for(unsigned int i = 0 ; i < size ; i++){
- pTreePath.push_back(lines[i]);
- }
+void RootInput::ReadOldStyleInputFile(NPL::InputParser& parser) {
+ vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("TTreeName");
+ pTreeName = blocks[0]->GetLines()[0];
+ vector<NPL::InputBlock*> files = parser.GetAllBlocksWithToken("RootFileName");
+ if (files.size() > 0) {
+ vector<string> lines = files[0]->GetLines();
+ unsigned int size = lines.size();
+ for (unsigned int i = 0; i < size; i++) {
+ pTreePath.push_back(lines[i]);
}
+ }
}
////////////////////////////////////////////////////////////////////////////////
-void RootInput::ReadInputFile(NPL::InputParser& parser){
- vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("Tree");
- pTreeName = blocks[0]->GetMainValue();
- std::vector<std::string> lines=blocks[0]->GetLines();
- unsigned int size = lines.size();
- for(unsigned int i = 0 ; i < size ; i++){
- if(lines[i].find(".root")!=string::npos)
- pTreePath.push_back(lines[i]);
- else if(lines[i].find(".tree")!=string::npos)
- ReadTreeFile(lines[i]);
- }
+void RootInput::ReadInputFile(NPL::InputParser& parser) {
+ vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("Tree");
+ pTreeName = blocks[0]->GetMainValue();
+ std::vector<std::string> lines = blocks[0]->GetLines();
+ unsigned int size = lines.size();
+ for (unsigned int i = 0; i < size; i++) {
+ if (lines[i].find(".root") != string::npos)
+ pTreePath.push_back(lines[i]);
+ else if (lines[i].find(".tree") != string::npos)
+ ReadTreeFile(lines[i]);
+ }
- vector<NPL::InputBlock*> friends = parser.GetAllBlocksWithToken("Friend");
- unsigned int sizeF = friends.size();
- for(unsigned int i = 0 ; i < sizeF ; i++){
- pFriendsPath.insert(pair< string,vector<string> > (friends[i]->GetMainValue(),friends[i]->GetLines()));
- }
+ vector<NPL::InputBlock*> friends = parser.GetAllBlocksWithToken("Friend");
+ unsigned int sizeF = friends.size();
+ for (unsigned int i = 0; i < sizeF; i++) {
+ pFriendsPath.insert(pair<string, vector<string>>(friends[i]->GetMainValue(), friends[i]->GetLines()));
+ }
}
////////////////////////////////////////////////////////////////////////////////
-void RootInput::ReadTreeFile(std::string path){
+void RootInput::ReadTreeFile(std::string path) {
ifstream tree(path.c_str());
- path=path.substr(0,path.rfind("/")+1);
+ path = path.substr(0, path.rfind("/") + 1);
std::string buffer;
- bool first=true;
- unsigned int count = 0 ;
- while(tree>>buffer){
- if(first){
- pTreePath.push_back(path+buffer);
+ bool first = true;
+ unsigned int count = 0;
+ while (tree >> buffer) {
+ if (first) {
+ pTreePath.push_back(path + buffer);
count++;
- first=false;
+ first = false;
}
- else{
- pFriendsTreePath[count++].push_back(path+buffer);
+ else {
+ pFriendsTreePath[count++].push_back(path + buffer);
}
}
}
////////////////////////////////////////////////////////////////////////////////
-RootInput::RootInput(std::string configFileName){
+RootInput::RootInput(std::string configFileName) {
NumberOfFriend = 0;
- pRootFile = NULL;
- std::string lastfile= NPOptionManager::getInstance()->GetLastFile();
- if(lastfile!="VOID"){
+ pRootFile = NULL;
+ std::string lastfile = NPOptionManager::getInstance()->GetLastFile();
+ if (lastfile != "VOID") {
configFileName = lastfile;
}
@@ -120,25 +120,25 @@ RootInput::RootInput(std::string configFileName){
// Old style file
vector<NPL::InputBlock*> old_blocks = parser.GetAllBlocksWithToken("TTreeName");
- if(old_blocks.size()==1){
+ if (old_blocks.size() == 1) {
ReadOldStyleInputFile(parser);
}
// New style file
vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("Tree");
- if(blocks.size()==1){
+ if (blocks.size() == 1) {
ReadInputFile(parser);
}
// If the tree come from a simulation, the InteractionCoordinates
// and InitialConditions lib are loaded
- if(pTreeName=="SimulatedTree"){
- std::string path = getenv("NPTOOL");
- path+="/NPLib/lib/";
- std::string libName="libNPInteractionCoordinates"+NPOptionManager::getInstance()->GetSharedLibExtension();
- libName=path+libName;
- dlopen(libName.c_str(),RTLD_NOW);
- libName="libNPInitialConditions"+NPOptionManager::getInstance()->GetSharedLibExtension();
- libName=path+libName;
- dlopen(libName.c_str(),RTLD_NOW);
+ if (pTreeName == "SimulatedTree") {
+ std::string path = getenv("NPTOOL");
+ path += "/NPLib/lib/";
+ std::string libName = "libNPInteractionCoordinates" + NPOptionManager::getInstance()->GetSharedLibExtension();
+ libName = path + libName;
+ dlopen(libName.c_str(), RTLD_NOW);
+ libName = "libNPInitialConditions" + NPOptionManager::getInstance()->GetSharedLibExtension();
+ libName = path + libName;
+ dlopen(libName.c_str(), RTLD_NOW);
}
// Initialise the chain
@@ -147,106 +147,111 @@ RootInput::RootInput(std::string configFileName){
// Add all the files
unsigned int size = pTreePath.size();
std::string firstfile;
- for(unsigned int i = 0 ; i < size ; i++){
+ for (unsigned int i = 0; i < size; i++) {
cout << " - Adding file : " << pTreePath[i].c_str() << endl;
pRootChain->Add(pTreePath[i].c_str());
// Test if the file is a regex or a single file
double counts;
std::string command = "ls " + pTreePath[i] + " > .ls_return";
- counts= system(command.c_str());
+ counts = system(command.c_str());
std::ifstream return_ls(".ls_return");
std::string files;
- while(return_ls >> files){
- if(counts == 0)
+ while (return_ls >> files) {
+ if (counts == 0)
firstfile = files;
counts++;
}
}
-
+
// Case of user made Friends
- for(auto it = pFriendsPath.begin(); it!=pFriendsPath.end() ; it++){
+ for (auto it = pFriendsPath.begin(); it != pFriendsPath.end(); it++) {
TChain* chain = new TChain(it->first.c_str());
cout << " - Adding friend : " << endl;
- for(auto itp = it->second.begin() ; itp!=it->second.end() ; itp++){
- cout << " - " << (*itp).c_str() << endl;
- chain->Add((*itp).c_str());
+ for (auto itp = it->second.begin(); itp != it->second.end(); itp++) {
+ cout << " - " << (*itp).c_str() << endl;
+ chain->Add((*itp).c_str());
}
pRootChain->AddFriend(chain);
}
// Case of tree file
- for(auto it = pFriendsTreePath.begin(); it!=pFriendsTreePath.end() ; it++){
+ for (auto it = pFriendsTreePath.begin(); it != pFriendsTreePath.end(); it++) {
TChain* chain = new TChain(pTreeName.c_str());
cout << " - Adding friend : " << endl;
- for(auto itp = it->second.begin() ; itp!=it->second.end() ; itp++){
+ for (auto itp = it->second.begin(); itp != it->second.end(); itp++) {
cout << " - " << (*itp).c_str() << endl;
chain->Add((*itp).c_str());
- }
+ }
pRootChain->AddFriend(chain);
}
-
- if (!pRootFile)
+
+ if (!pRootFile)
pRootFile = new TFile(firstfile.c_str());
- if( pRootChain->GetEntries() ==0){
- std::cout << "\033[1;31m**** ERROR: No entries to analyse ****\033[0m" << std::endl;
+ if (pRootChain->GetEntries() == 0) {
+ std::cout << "\033[1;31m**** ERROR: No entries to analyse ****\033[0m" << std::endl;
exit(1);
}
- else
- std::cout << "\033[1;32mROOTInput: " << pRootChain->GetEntries() << " entries loaded in the input chain\033[0m" << std::endl ;
+ else
+ std::cout << "\033[1;32mROOTInput: " << pRootChain->GetEntries() << " entries loaded in the input chain\033[0m"
+ << std::endl;
}
////////////////////////////////////////////////////////////////////////////////
-std::string RootInput::DumpAsciiFile(const char* type, const char* folder){
- std::string name="fail";
+std::string RootInput::DumpAsciiFile(const char* type, const char* folder) {
+ std::string name = "fail";
std::string sfolder = folder;
// create folder if not existing
// first test if folder exist
struct stat dirInfo;
int res = stat(folder, &dirInfo);
- if (res != 0) { // if folder does not exist, create it
+ if (res != 0) { // if folder does not exist, create it
std::string cmd = "mkdir " + sfolder;
- if (system(cmd.c_str()) != 0) std::cout << "RootInput::DumpAsciiFile() problem creating directory" << std::endl;
+ if (system(cmd.c_str()) != 0)
+ std::cout << "RootInput::DumpAsciiFile() problem creating directory" << std::endl;
}
std::string stype = type;
if (stype == "EventGenerator") {
- TAsciiFile *aFile = (TAsciiFile*)pRootFile->Get(stype.c_str());
+ TAsciiFile* aFile = (TAsciiFile*)pRootFile->Get(stype.c_str());
- if(aFile)
- {
+ if (aFile) {
// build file name
std::string title = aFile->GetTitle();
size_t pos = title.rfind("/");
- if (pos != std::string::npos) name = sfolder + title.substr(pos);
- else name = sfolder + "/" + title;
+ if (pos != std::string::npos)
+ name = sfolder + title.substr(pos);
+ else
+ name = sfolder + "/" + title;
aFile->WriteToFile(name.c_str());
}
}
else if (stype == "DetectorConfiguration") {
- TAsciiFile *aFile = (TAsciiFile*)pRootFile->Get(stype.c_str());
- if(aFile)
- {
+ TAsciiFile* aFile = (TAsciiFile*)pRootFile->Get(stype.c_str());
+ if (aFile) {
// build file name
std::string title = aFile->GetTitle();
size_t pos = title.rfind("/");
- if (pos != std::string::npos) name = sfolder + title.substr(pos);
- else name = sfolder + "/" + title;
+ if (pos != std::string::npos)
+ name = sfolder + title.substr(pos);
+ else
+ name = sfolder + "/" + title;
aFile->WriteToFile(name.c_str());
}
}
else if (stype == "Calibration") {
- TAsciiFile *aFile = (TAsciiFile*)pRootFile->Get(stype.c_str());
- if(aFile)
- {
+ TAsciiFile* aFile = (TAsciiFile*)pRootFile->Get(stype.c_str());
+ if (aFile) {
// build file name
std::string title = aFile->GetTitle();
size_t pos = title.rfind("/");
- if (pos != std::string::npos) name = sfolder + title.substr(pos);
- else name = sfolder + "/" + title;
+ if (pos != std::string::npos)
+ name = sfolder + title.substr(pos);
+ else
+ name = sfolder + "/" + title;
aFile->WriteToFile(name.c_str());
}
}
@@ -261,12 +266,13 @@ std::string RootInput::DumpAsciiFile(const char* type, const char* folder){
}
////////////////////////////////////////////////////////////////////////////////
-RootInput::~RootInput(){
+RootInput::~RootInput() {
// delete default directory ./.tmp
struct stat dirInfo;
int res = stat("./.tmp", &dirInfo);
- if (res == 0) { // if does exist, delete it
- if (system("rm -rf ./.tmp") != 0) std::cout << "RootInput::~RootInput() problem deleting ./.tmp directory" << std::endl;
+ if (res == 0) { // if does exist, delete it
+ if (system("rm -rf ./.tmp") != 0)
+ std::cout << "RootInput::~RootInput() problem deleting ./.tmp directory" << std::endl;
}
std::cout << std::endl << "Root Input summary" << std::endl;
std::cout << " - Number of bites read: " << pRootFile->GetBytesRead() << std::endl;
diff --git a/NPLib/Detectors/MDM/CMakeLists.txt b/NPLib/Detectors/MDM/CMakeLists.txt
index a3ad822b2..6d3113af0 100644
--- a/NPLib/Detectors/MDM/CMakeLists.txt
+++ b/NPLib/Detectors/MDM/CMakeLists.txt
@@ -6,7 +6,7 @@ add_custom_command(OUTPUT TMDMDataDict.cxx COMMAND ${CMAKE_BINARY_DIR}/scripts/b
##
## Check if MINUIT2 is installed along with ROOT
-find_library(libMinuit2_FOUND NAMES Minuit2 HINTS "${ROOTSYS}/lib")
+find_library(libMinuit2_FOUND NAMES Minuit2 HINTS "${ROOTSYS}/lib" "/usr/lib64/root")
if(libMinuit2_FOUND)
message(STATUS "Minuit2 support enabled for MDM.")
add_definitions(-DHAVE_MINUIT2)
diff --git a/NPLib/Detectors/Minos/CMakeLists.txt b/NPLib/Detectors/Minos/CMakeLists.txt
index d81fe157d..6f6affcaa 100644
--- a/NPLib/Detectors/Minos/CMakeLists.txt
+++ b/NPLib/Detectors/Minos/CMakeLists.txt
@@ -10,7 +10,7 @@ add_custom_command(OUTPUT TMinosClustDict.cxx COMMAND ../../scripts/build_dict.s
add_custom_command(OUTPUT TMinosResultDict.cxx COMMAND ../../scripts/build_dict.sh TMinosResult.h TMinosResultDict.cxx TMinosResult.rootmap libNPMinos.dylib DEPENDS TMinosResult.h)
## Check if MINUIT is installed along with ROOT
-find_library(libMinuit_FOUND NAMES Minuit2 HINTS "$ENV{ROOTSYS}/lib")
+find_library(libMinuit_FOUND NAMES Minuit2 HINTS "$ENV{ROOTSYS}/lib" "/usr/lib64/root")
if(libMinuit_FOUND)
message(STATUS "Minuit support enabled for Minos.")
add_definitions(-DHAVE_MINUIT)
diff --git a/NPLib/Detectors/Minos/TMinosPhysics.cxx b/NPLib/Detectors/Minos/TMinosPhysics.cxx
index cc2e729e7..04ccd7800 100644
--- a/NPLib/Detectors/Minos/TMinosPhysics.cxx
+++ b/NPLib/Detectors/Minos/TMinosPhysics.cxx
@@ -10,33 +10,33 @@
* *
* Creation Date : 2014/11/24 *
* Last update : 2019/09 implemeted in NPTool by Cyril Lenain *
- * lenain@lpccaen.in2p3.fr *
+ * lenain@lpccaen.in2p3.fr *
*---------------------------------------------------------------------------*
* Decription: *
* This class hold Minos Treated data *
* *
*---------------------------------------------------------------------------*
* Comment: *
- * *
+ * *
* *
*****************************************************************************/
#include "TMinosPhysics.h"
// STL
-#include <sstream>
-#include <iostream>
#include <cmath>
-#include <stdlib.h>
+#include <iostream>
#include <limits>
+#include <sstream>
+#include <stdlib.h>
using namespace std;
// NPL
-#include "RootInput.h"
-#include "RootOutput.h"
#include "NPDetectorFactory.h"
#include "NPOptionManager.h"
#include "NPSystemOfUnits.h"
+#include "RootInput.h"
+#include "RootOutput.h"
// ROOT
#include "TChain.h"
@@ -44,100 +44,94 @@ using namespace std;
TMinosPhysics* current_phy = 0;
///////////////////////////////////////////////////////////////////////////
-TMinosPhysics::TMinosPhysics(){
- m_EventData=new TMinosData;
- m_EventPhysics=this;
- m_E_RAW_Threshold=0; // adc channels
- m_E_Threshold=0; // MeV
- m_ransac.SetParameters(100,// max 100 iteration
- 5, // a point belong to a track if it is within 5 mm
- 40, // max distance allowed between a pair of point to consider a track
- 10);// minimum point to form a valid cluster
+TMinosPhysics::TMinosPhysics() {
+ m_EventData = new TMinosData;
+ m_EventPhysics = this;
+ m_E_RAW_Threshold = 0; // adc channels
+ m_E_Threshold = 0; // MeV
+ m_ransac.SetParameters(100, // max 100 iteration
+ 5, // a point belong to a track if it is within 5 mm
+ 40, // max distance allowed between a pair of point to consider a track
+ 10); // minimum point to form a valid cluster
}
///////////////////////////////////////////////////////////////////////////
-void TMinosPhysics::BuildSimplePhysicalEvent() {
- BuildPhysicalEvent();
-}
+void TMinosPhysics::BuildSimplePhysicalEvent() { BuildPhysicalEvent(); }
///////////////////////////////////////////////////////////////////////////
void TMinosPhysics::BuildPhysicalEvent() {
PreTreat();
- vector<NPL::LinearCluster3D> clusters = m_ransac.TreatEvent(X_Pad,Y_Pad,Z_Pad);
+ vector<NPL::LinearCluster3D> clusters = m_ransac.TreatEvent(X_Pad, Y_Pad, Z_Pad);
unsigned int sizeC = clusters.size();
- for(unsigned int i = 0 ; i < sizeC ; i++){
- // Extract the Direction of the track
+ for (unsigned int i = 0; i < sizeC; i++) {
+ // Extract the Direction of the track
clusters[i].LinearFit();
Tracks_P0.push_back(clusters[i].GetP0());
Tracks_Dir.push_back(clusters[i].GetDir());
}
- if(sizeC==2){
- static TVector3 Vertex,delta;
- Delta_Vertex =
- MinimumDistanceTwoLines(Tracks_P0[0],Tracks_P0[0]+Tracks_Dir[0],
- Tracks_P0[1],Tracks_P0[1]+Tracks_Dir[1],
- Vertex, delta);
- X_Vertex= Vertex.X();
- Y_Vertex= Vertex.Y();
- Z_Vertex= Vertex.Z();
- Theta_12=180.*Tracks_Dir[0].Angle(Tracks_Dir[1])/(M_PI);
+ if (sizeC == 2) {
+ static TVector3 Vertex, delta;
+ Delta_Vertex = MinimumDistanceTwoLines(Tracks_P0[0], Tracks_P0[0] + Tracks_Dir[0], Tracks_P0[1],
+ Tracks_P0[1] + Tracks_Dir[1], Vertex, delta);
+ X_Vertex = Vertex.X();
+ Y_Vertex = Vertex.Y();
+ Z_Vertex = Vertex.Z();
+ Theta_12 = 180. * Tracks_Dir[0].Angle(Tracks_Dir[1]) / (M_PI);
}
}
///////////////////////////////////////////////////////////////////////////
void TMinosPhysics::PreTreat() {
// apply thresholds and calibration
- static unsigned int sizePad,sizeQ ;
+ static unsigned int sizePad, sizeQ;
static unsigned short PadNumber;
- static double Q,T;
- static auto cal= CalibrationManager::getInstance();
+ static double Q, T;
+ static auto cal = CalibrationManager::getInstance();
static string cal_v, cal_o;
sizePad = m_EventData->GetPadMult();
- if(sizePad>20){
- for(unsigned int i = 0 ; i < sizePad ; i++){
+ if (sizePad > 20) {
+ for (unsigned int i = 0; i < sizePad; i++) {
vector<unsigned short>* Charge = m_EventData->GetChargePtr(i);
- if(Charge->size()>40 ){
- vector<unsigned short>* Time = m_EventData->GetTimePtr(i);
- m_utility.Calibrate(Time,Charge,i,T,Q);
+ if (Charge->size() > 40) {
+ vector<unsigned short>* Time = m_EventData->GetTimePtr(i);
+ m_utility.Calibrate(Time, Charge, i, T, Q);
- if(T>0&&Q<77000){
+ if (T > 0 && Q < 77000) {
PadNumber = m_EventData->GetPadNumber(i);
double x_mm = m_X[PadNumber];
double y_mm = m_Y[PadNumber];
- unsigned int ring = round((sqrt(x_mm*x_mm + y_mm*y_mm)-44.15)/2.1);
- cal_v="Minos/R"+NPL::itoa(ring)+"_VDRIFT";
- cal_o="Minos/R"+NPL::itoa(ring)+"_OFFSET";
- double calV= cal->GetValue(cal_v,0);
- double calO=cal->GetValue(cal_o,0);
-
- double z_mm = (T*m_TimeBin+calO/*-0.7-0.04*ring*/)*calV;
- //cout << T*m_TimeBin+calO << endl;
- TVector3 Pos=TVector3(x_mm+m_Position.X(),y_mm+m_Position.Y(),z_mm+m_Position.Z());
- Pos.RotateZ(m_ZRotation);
+ unsigned int ring = round((sqrt(x_mm * x_mm + y_mm * y_mm) - 44.15) / 2.1);
+ cal_v = "Minos/R" + NPL::itoa(ring) + "_VDRIFT";
+ cal_o = "Minos/R" + NPL::itoa(ring) + "_OFFSET";
+ double calV = cal->GetValue(cal_v, 0);
+ double calO = cal->GetValue(cal_o, 0);
+
+ double z_mm = (T * m_TimeBin + calO /*-0.7-0.04*ring*/) * calV;
+ // cout << T*m_TimeBin+calO << endl;
+ TVector3 Pos = TVector3(x_mm + m_Position.X(), y_mm + m_Position.Y(), z_mm + m_Position.Z());
+ Pos.RotateZ(m_ZRotation);
// Calibrate the Pad:
-
+
X_Pad.push_back(Pos.X());
Y_Pad.push_back(Pos.Y());
- Z_Pad.push_back(Pos.Z());
+ Z_Pad.push_back(Pos.Z());
Ring_Pad.push_back(ring);
- Q_Pad.push_back(Q);
- T_Pad.push_back(T);
+ Q_Pad.push_back(Q);
+ T_Pad.push_back(T);
}
}
- else{
+ else {
/* X_Pad.push_back(-1);
Y_Pad.push_back(-1);
Z_Pad.push_back((-1);
Ring_Pad.push_back(-1);
- Z_Pad.push_back(-1);
- Q_Pad.push_back(-1);
- T_Pad.push_back(-1);*/
+ Z_Pad.push_back(-1);
+ Q_Pad.push_back(-1);
+ T_Pad.push_back(-1);*/
}
-
}
-
}
}
///////////////////////////////////////////////////////////////////////////
@@ -162,33 +156,33 @@ void TMinosPhysics::ReadAnalysisConfig() {
asciiConfig->Append(FileName.c_str());
asciiConfig->AppendLine("");
// read analysis config file
- string LineBuffer,DataBuffer,whatToDo;
+ string LineBuffer, DataBuffer, whatToDo;
while (!AnalysisConfigFile.eof()) {
// Pick-up next line
getline(AnalysisConfigFile, LineBuffer);
// search for "header"
string name = "ConfigMinos";
- if (LineBuffer.compare(0, name.length(), name) == 0)
+ if (LineBuffer.compare(0, name.length(), name) == 0)
ReadingStatus = true;
// loop on tokens and data
- while (ReadingStatus ) {
- whatToDo="";
+ while (ReadingStatus) {
+ whatToDo = "";
AnalysisConfigFile >> whatToDo;
// Search for comment symbol (%)
if (whatToDo.compare(0, 1, "%") == 0) {
- AnalysisConfigFile.ignore(numeric_limits<streamsize>::max(), '\n' );
+ AnalysisConfigFile.ignore(numeric_limits<streamsize>::max(), '\n');
}
- else if (whatToDo=="E_RAW_THRESHOLD") {
+ else if (whatToDo == "E_RAW_THRESHOLD") {
AnalysisConfigFile >> DataBuffer;
m_E_RAW_Threshold = atof(DataBuffer.c_str());
cout << whatToDo << " " << m_E_RAW_Threshold << endl;
}
- else if (whatToDo=="E_THRESHOLD") {
+ else if (whatToDo == "E_THRESHOLD") {
AnalysisConfigFile >> DataBuffer;
m_E_Threshold = atof(DataBuffer.c_str());
cout << whatToDo << " " << m_E_Threshold << endl;
@@ -199,14 +193,13 @@ void TMinosPhysics::ReadAnalysisConfig() {
}
}
}
-
}
///////////////////////////////////////////////////////////////////////////
void TMinosPhysics::Clear() {
- Ring_Pad.clear();
- X_Pad.clear();
- Y_Pad.clear();
+ Ring_Pad.clear();
+ X_Pad.clear();
+ Y_Pad.clear();
Z_Pad.clear();
Q_Pad.clear();
T_Pad.clear();
@@ -216,91 +209,86 @@ void TMinosPhysics::Clear() {
Tracks_Dir.clear();
// Vertex information
- X_Vertex=-1000;
- Y_Vertex=-1000;
- Z_Vertex=-1000;
- Theta_12=-1000;
- Delta_Vertex=-1000;
+ X_Vertex = -1000;
+ Y_Vertex = -1000;
+ Z_Vertex = -1000;
+ Theta_12 = -1000;
+ Delta_Vertex = -1000;
}
///////////////////////////////////////////////////////////////////////////
void TMinosPhysics::ReadConfiguration(NPL::InputParser parser) {
vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("Minos");
- if(NPOptionManager::getInstance()->GetVerboseLevel())
- cout << "//// " << blocks.size() << " detector(s) found " << endl;
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << "//// " << blocks.size() << " detector(s) found " << endl;
- vector<string> token= {"XML","TimeBin","ShapingTime","Baseline","Sampling","Position","ZRotation"};
+ vector<string> token = {"XML", "TimeBin", "ShapingTime", "Baseline", "Sampling", "Position", "ZRotation"};
+
+ for (unsigned int i = 0; i < blocks.size(); i++) {
- for(unsigned int i = 0 ; i < blocks.size() ; i++){
-
if (blocks[i]->HasTokenList(token)) {
- cout << endl << "//// MINOS (" << i+1 << ")" << endl;
- m_ShapingTime = blocks[i]->GetDouble("ShapingTime","ns")/NPUNITS::ns;
- m_TimeBin = blocks[i]->GetDouble("TimeBin","ns")/NPUNITS::ns;
- m_Sampling= blocks[i]->GetInt("Sampling");
- m_Baseline= blocks[i]->GetInt("BaseLine");
- m_utility.SetParameters(m_TimeBin,m_ShapingTime,m_Baseline,m_Sampling);
- m_Position = blocks[i]->GetTVector3("Position","mm");
- m_ZRotation= blocks[i]->GetDouble("ZRotation","deg");
+ cout << endl << "//// MINOS (" << i + 1 << ")" << endl;
+ m_ShapingTime = blocks[i]->GetDouble("ShapingTime", "ns") / NPUNITS::ns;
+ m_TimeBin = blocks[i]->GetDouble("TimeBin", "ns") / NPUNITS::ns;
+ m_Sampling = blocks[i]->GetInt("Sampling");
+ m_Baseline = blocks[i]->GetInt("BaseLine");
+ m_utility.SetParameters(m_TimeBin, m_ShapingTime, m_Baseline, m_Sampling);
+ m_Position = blocks[i]->GetTVector3("Position", "mm");
+ m_ZRotation = blocks[i]->GetDouble("ZRotation", "deg");
string xmlpath = blocks[i]->GetString("XML");
NPL::XmlParser xml;
xml.LoadFile(xmlpath);
ReadXML(xml);
-
}
else {
- cout << "ERROR: Missing token for Minos, check your input file"<< endl;
- cout << "Required token: ";
- for(unsigned int i = 0 ; i < token.size() ; i++)
- cout << token[i] << " " ;
+ cout << "ERROR: Missing token for Minos, check your input file" << endl;
+ cout << "Required token: ";
+ for (unsigned int i = 0; i < token.size(); i++)
+ cout << token[i] << " ";
cout << endl;
exit(1);
}
-
}
}
///////////////////////////////////////////////////////////////////////////
-void TMinosPhysics::ReadXML(NPL::XmlParser& xml){
- std::vector<NPL::XML::block*> b = xml.GetAllBlocksWithName("MINOS");
+void TMinosPhysics::ReadXML(NPL::XmlParser& xml) {
+ std::vector<NPL::XML::block*> b = xml.GetAllBlocksWithName("MINOS");
unsigned int size = b.size();
- unsigned int valid_pad=0;
- for(unsigned int i = 0 ; i < size ; i++){
+ unsigned int valid_pad = 0;
+ for (unsigned int i = 0; i < size; i++) {
- unsigned short ID = b[i]->AsInt("ID");
- m_X[ID] = b[i]->AsDouble("X");
- m_Y[ID] = b[i]->AsDouble("Y");
+ unsigned short ID = b[i]->AsInt("ID");
+ m_X[ID] = b[i]->AsDouble("X");
+ m_Y[ID] = b[i]->AsDouble("Y");
valid_pad++;
- if((m_X[ID]==0&&m_Y[ID]==0)||(m_X[ID]==-1&&m_Y[ID]==-1))
+ if ((m_X[ID] == 0 && m_Y[ID] == 0) || (m_X[ID] == -1 && m_Y[ID] == -1))
valid_pad--;
- m_Period[ID] = b[i]->AsDouble("TPERIOD");
- m_Offset[ID] = b[i]->AsDouble("TOFFSET");
- m_Gain[ID] = b[i]->AsDouble("GAIN");
+ m_Period[ID] = b[i]->AsDouble("TPERIOD");
+ m_Offset[ID] = b[i]->AsDouble("TOFFSET");
+ m_Gain[ID] = b[i]->AsDouble("GAIN");
}
- if(NPOptionManager::getInstance()->GetVerboseLevel())
- cout << "//// " << m_X.size() << " pads found with " << valid_pad << " valid pads" << endl;
-
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << "//// " << m_X.size() << " pads found with " << valid_pad << " valid pads" << endl;
}
-
-
///////////////////////////////////////////////////////////////////////////
void TMinosPhysics::AddParameterToCalibrationManager() {
CalibrationManager* Cal = CalibrationManager::getInstance();
unsigned int NbrRing = 18;
- vector<double> standardV={0.03475};
- vector<double> standardO={0};
- for(unsigned int i = 0 ; i < NbrRing ; i++){
- Cal->AddParameter("Minos", "R"+NPL::itoa(i+1)+"_VDRIFT","Minos_R"+NPL::itoa(i+1)+"_VDRIFT",standardV);
- Cal->AddParameter("Minos", "R"+NPL::itoa(i+1)+"_OFFSET","Minos_R"+NPL::itoa(i+1)+"_OFFSET",standardO);
+ vector<double> standardV = {0.03475};
+ vector<double> standardO = {0};
+ for (unsigned int i = 0; i < NbrRing; i++) {
+ Cal->AddParameter("Minos", "R" + NPL::itoa(i + 1) + "_VDRIFT", "Minos_R" + NPL::itoa(i + 1) + "_VDRIFT", standardV);
+ Cal->AddParameter("Minos", "R" + NPL::itoa(i + 1) + "_OFFSET", "Minos_R" + NPL::itoa(i + 1) + "_OFFSET", standardO);
}
}
///////////////////////////////////////////////////////////////////////////
void TMinosPhysics::InitializeRootInputRaw() {
TChain* inputChain = RootInput::getInstance()->GetChain();
- inputChain->SetBranchStatus("Minos", true );
- inputChain->SetBranchAddress("Minos", &m_EventData );
+ inputChain->SetBranchStatus("Minos", true);
+ inputChain->SetBranchAddress("Minos", &m_EventData);
}
///////////////////////////////////////////////////////////////////////////
@@ -318,21 +306,19 @@ void TMinosPhysics::InitializeRootOutput() {
////////////////////////////////////////////////////////////////////////////////
// Construct Method to be pass to the DetectorFactory //
////////////////////////////////////////////////////////////////////////////////
-NPL::VDetector* TMinosPhysics::Construct() {
- return (NPL::VDetector*) new TMinosPhysics();
-}
+NPL::VDetector* TMinosPhysics::Construct() { return (NPL::VDetector*)new TMinosPhysics(); }
////////////////////////////////////////////////////////////////////////////////
// Registering the construct method to the factory //
////////////////////////////////////////////////////////////////////////////////
-extern "C"{
- class proxy_Minos{
- public:
- proxy_Minos(){
- NPL::DetectorFactory::getInstance()->AddToken("Minos","Minos");
- NPL::DetectorFactory::getInstance()->AddDetector("Minos",TMinosPhysics::Construct);
- }
- };
+extern "C" {
+class proxy_Minos {
+ public:
+ proxy_Minos() {
+ NPL::DetectorFactory::getInstance()->AddToken("Minos", "Minos");
+ NPL::DetectorFactory::getInstance()->AddDetector("Minos", TMinosPhysics::Construct);
+ }
+};
- proxy_Minos p_Minos;
+proxy_Minos p_Minos;
}
diff --git a/NPLib/TrackReconstruction/CMakeLists.txt b/NPLib/TrackReconstruction/CMakeLists.txt
index 341a03e34..f87460b71 100644
--- a/NPLib/TrackReconstruction/CMakeLists.txt
+++ b/NPLib/TrackReconstruction/CMakeLists.txt
@@ -7,7 +7,7 @@ add_custom_command(OUTPUT NPClusterDict.cxx COMMAND ${CMAKE_BINARY_DIR}/scripts/
add_custom_command(OUTPUT TrackingDict.cxx COMMAND ${CMAKE_BINARY_DIR}/scripts/build_dict.sh Tracking.h TrackingDict.cxx Tracking.rootmap libNPTrackReconstruction.so NPTrackReconstructionLinkDef.h DEPENDS Tracking.h)
## Check if MINUIT2 is installed along with ROOT
-find_library(libMinuit2_FOUND NAMES Minuit2 HINTS "${ROOTSYS}/lib")
+find_library(libMinuit2_FOUND NAMES Minuit2 HINTS "${ROOTSYS}/lib" "/usr/lib64/root")
if(libMinuit2_FOUND)
message(STATUS "Minuit2 support enabled for TrackReconstruction.")
add_definitions(-DHAVE_MINUIT2)
diff --git a/NPLib/TrackReconstruction/NPLinearRansac3D.cxx b/NPLib/TrackReconstruction/NPLinearRansac3D.cxx
index 767e43ce0..8ef0839a4 100644
--- a/NPLib/TrackReconstruction/NPLinearRansac3D.cxx
+++ b/NPLib/TrackReconstruction/NPLinearRansac3D.cxx
@@ -16,139 +16,143 @@
* Find linear tracks using the ransac method *
*****************************************************************************/
-#include"NPLinearRansac3D.h"
+#include "NPLinearRansac3D.h"
using namespace std;
using namespace NPL;
////////////////////////////////////////////////////////////////////////////////
-void LinearCluster3D::LinearFit(){
+void LinearCluster3D::LinearFit() {
- unsigned int sizeI = m_index.size();
+ unsigned int sizeI = m_index.size();
double min_distance = 1e20;
- unsigned int min_i = 0 ; unsigned int min_j = 0;
-
- double meanDX=0;
- double meanDY=0;
- double meanDZ=0;
- double totalW=0;
- double meanX=0;
- double meanY=0;
- double meanZ=0;
+ unsigned int min_i = 0;
+ unsigned int min_j = 0;
+
+ double meanDX = 0;
+ double meanDY = 0;
+ double meanDZ = 0;
+ double totalW = 0;
+ double meanX = 0;
+ double meanY = 0;
+ double meanZ = 0;
double w;
- for(unsigned int i = 0 ; i < sizeI ; i++){
- TVector3 Pi((*m_X)[m_index[i]],(*m_Y)[m_index[i]],(*m_Z)[m_index[i]]);
- for(unsigned int j = i+1 ; j < sizeI ; j++){
- double dij=0;
- TVector3 Pj((*m_X)[m_index[j]],(*m_Y)[m_index[j]],(*m_Z)[m_index[j]]);
+ for (unsigned int i = 0; i < sizeI; i++) {
+ TVector3 Pi((*m_X)[m_index[i]], (*m_Y)[m_index[i]], (*m_Z)[m_index[i]]);
+ for (unsigned int j = i + 1; j < sizeI; j++) {
+ double dij = 0;
+ TVector3 Pj((*m_X)[m_index[j]], (*m_Y)[m_index[j]], (*m_Z)[m_index[j]]);
// compute the average distance of all point to this line
- for(unsigned int p = 0 ; p < sizeI ; p++){
- TVector3 Pp((*m_X)[m_index[p]],(*m_Y)[m_index[p]],(*m_Z)[m_index[p]]);
- dij+=MinimumDistancePointLine(Pi,Pj,Pp);
+ for (unsigned int p = 0; p < sizeI; p++) {
+ TVector3 Pp((*m_X)[m_index[p]], (*m_Y)[m_index[p]], (*m_Z)[m_index[p]]);
+ dij += MinimumDistancePointLine(Pi, Pj, Pp);
}
-
- w = 1./dij;
- totalW+=w;
-
- meanX+=w*((*m_X)[m_index[i]]);
- meanY+=w*((*m_Y)[m_index[i]]);
- meanZ+=w*((*m_Z)[m_index[i]]);
-
- meanDX+=w*((*m_X)[m_index[i]]-(*m_X)[m_index[j]]);
- meanDY+=w*((*m_Y)[m_index[i]]-(*m_Y)[m_index[j]]);
- meanDZ+=w*((*m_Z)[m_index[i]]-(*m_Z)[m_index[j]]);
- }
- }
- meanDX/=totalW;
- meanDY/=totalW;
- meanDZ/=totalW;
- meanX/=totalW;
- meanY/=totalW;
- meanZ/=totalW;
+ w = 1. / dij;
+ totalW += w;
- m_P0 = TVector3(meanX,meanY,meanZ);
- m_Dir = (TVector3(meanDX,meanDY,meanDZ)).Unit();
+ meanX += w * ((*m_X)[m_index[i]]);
+ meanY += w * ((*m_Y)[m_index[i]]);
+ meanZ += w * ((*m_Z)[m_index[i]]);
+ meanDX += w * ((*m_X)[m_index[i]] - (*m_X)[m_index[j]]);
+ meanDY += w * ((*m_Y)[m_index[i]] - (*m_Y)[m_index[j]]);
+ meanDZ += w * ((*m_Z)[m_index[i]] - (*m_Z)[m_index[j]]);
+ }
+ }
+
+ meanDX /= totalW;
+ meanDY /= totalW;
+ meanDZ /= totalW;
+ meanX /= totalW;
+ meanY /= totalW;
+ meanZ /= totalW;
+
+ m_P0 = TVector3(meanX, meanY, meanZ);
+ m_Dir = (TVector3(meanDX, meanDY, meanDZ)).Unit();
+ // Cluster dir is arbitrary. we choose to always have positive Z dir
+ if (m_Dir.Z() < 0)
+ m_Dir = -m_Dir;
};
////////////////////////////////////////////////////////////////////////////////
- vector<LinearCluster3D> LinearRansac3D::TreatEvent(vector<double>& X, vector<double>&Y, vector<double>&Z){
- cluster_id.clear();
- clusters.clear();
- unsigned int sizeX = X.size();
- cluster_id.resize(sizeX,0);
- m_cluster.clear();
- m_assigned.clear();
-
- if(sizeX<m_min_count)
- return clusters;
-
- unsigned int p1,p2;
- double d;
- TVector3 Vp1,Vp2,D,P;
- m_iteration_d=0;
- m_iteration=0;
-
- while(m_iteration++ < m_max_iteration && m_assigned.size()<sizeX){
- LinearCluster3D cluster(&X,&Y,&Z);
- // take 2 distant point randomly that has not been match before
- d=0 ; m_iteration_d=0;
- while(d<3*m_match_distance && m_iteration_d++<m_max_iteration ){
- p1 = rand.Integer(sizeX);
- p2 = rand.Integer(sizeX);
-
- Vp1.SetXYZ(X[p1],Y[p1],Z[p1]);
- Vp2.SetXYZ(X[p2],Y[p2],Z[p2]);
- D=Vp1-Vp2;
- d = D.Mag();
- if(d>m_max_distance)
- d=0;
- }
+vector<LinearCluster3D> LinearRansac3D::TreatEvent(vector<double>& X, vector<double>& Y, vector<double>& Z) {
+ cluster_id.clear();
+ clusters.clear();
+ unsigned int sizeX = X.size();
+ cluster_id.resize(sizeX, 0);
+ m_cluster.clear();
+ m_assigned.clear();
+
+ if (sizeX < m_min_count)
+ return clusters;
+
+ unsigned int p1, p2;
+ double d;
+ TVector3 Vp1, Vp2, D, P;
+ m_iteration_d = 0;
+ m_iteration = 0;
+
+ while (m_iteration++ < m_max_iteration && m_assigned.size() < sizeX) {
+ LinearCluster3D cluster(&X, &Y, &Z);
+ // take 2 distant point randomly that has not been match before
+ d = 0;
+ m_iteration_d = 0;
+ while (d < 3 * m_match_distance && m_iteration_d++ < m_max_iteration) {
+ p1 = rand.Integer(sizeX);
+ p2 = rand.Integer(sizeX);
+
+ Vp1.SetXYZ(X[p1], Y[p1], Z[p1]);
+ Vp2.SetXYZ(X[p2], Y[p2], Z[p2]);
+ D = Vp1 - Vp2;
+ d = D.Mag();
+ if (d > m_max_distance)
+ d = 0;
+ }
- // loop over all points
- for(unsigned int i = 0 ; i < sizeX ; i++){
- P.SetXYZ(X[i],Y[i],Z[i]);
- if(MinimumDistancePointLine(Vp1,Vp2,P) < m_match_distance){
- cluster.AddIndex(i);
- m_assigned.insert(i);
- }
- }
- // insert the newly formed cluster
- if(cluster.size()>m_min_count){
- m_cluster.insert(cluster);
- }
+ // loop over all points
+ for (unsigned int i = 0; i < sizeX; i++) {
+ P.SetXYZ(X[i], Y[i], Z[i]);
+ if (MinimumDistancePointLine(Vp1, Vp2, P) < m_match_distance) {
+ cluster.AddIndex(i);
+ m_assigned.insert(i);
+ }
+ }
+ // insert the newly formed cluster
+ if (cluster.size() > m_min_count) {
+ m_cluster.insert(cluster);
+ }
- }//while
-
- // loop over the cluster starting with the biggest
- unsigned int current_cluster=0;
- unsigned int index;
- for(auto it = m_cluster.begin() ; it!=m_cluster.end() ; ++it){
- current_cluster++;
- // cout << current_cluster << endl;
- unsigned int sizeC = (*it).size();
- unsigned int cluster_size = 0;
- for(unsigned int i = 0 ; i < sizeC ; i++){
- // Assigned cluster id to identified points
- unsigned int index = (*it).GetIndex(i);
- if(!cluster_id[index]){
- cluster_id[index]=current_cluster;
- cluster_size++;
- }
- }
- if(cluster_size<m_min_count){
- for(unsigned int i = 0 ; i < sizeC ; i++){
- unsigned int index = (*it).GetIndex(i);
- // remove the assigned point
- if(cluster_id[index]==current_cluster){
- cluster_id[index]=0;
- }
- }
- }
- else{
- clusters.push_back(*it);
+ } // while
+
+ // loop over the cluster starting with the biggest
+ unsigned int current_cluster = 0;
+ unsigned int index;
+ for (auto it = m_cluster.begin(); it != m_cluster.end(); ++it) {
+ current_cluster++;
+ // cout << current_cluster << endl;
+ unsigned int sizeC = (*it).size();
+ unsigned int cluster_size = 0;
+ for (unsigned int i = 0; i < sizeC; i++) {
+ // Assigned cluster id to identified points
+ unsigned int index = (*it).GetIndex(i);
+ if (!cluster_id[index]) {
+ cluster_id[index] = current_cluster;
+ cluster_size++;
+ }
+ }
+ if (cluster_size < m_min_count) {
+ for (unsigned int i = 0; i < sizeC; i++) {
+ unsigned int index = (*it).GetIndex(i);
+ // remove the assigned point
+ if (cluster_id[index] == current_cluster) {
+ cluster_id[index] = 0;
}
}
- return clusters;
}
+ else {
+ clusters.push_back(*it);
+ }
+ }
+ return clusters;
+}
diff --git a/NPSimulation/Detectors/SuperX3/SuperX3.cc b/NPSimulation/Detectors/SuperX3/SuperX3.cc
index a571962c0..8c0ab28f2 100644
--- a/NPSimulation/Detectors/SuperX3/SuperX3.cc
+++ b/NPSimulation/Detectors/SuperX3/SuperX3.cc
@@ -183,7 +183,7 @@ void SuperX3::ReadConfiguration(NPL::InputParser parser) {
if (blocks[i]->GetInt("VIS"))
m_non_sensitive_part_visiualisation = true;
- AddDetector(Theta, Phi, R, beta[0], beta[1], beta[2]);
+ AddDetector(R, Theta, Phi, beta[0], beta[1], beta[2]);
}
else {
@@ -237,18 +237,18 @@ void SuperX3::ConstructDetector(G4LogicalVolume* world) {
// w perpendicular to (u,v) plan and pointing ThirdStage
// Phi is angle between X axis and projection in (X,Y) plan
// Theta is angle between position vector and z axis
- G4double wX = m_R[i] * sin(Theta / rad) * cos(Phi / rad);
- G4double wY = m_R[i] * sin(Theta / rad) * sin(Phi / rad);
- G4double wZ = m_R[i] * cos(Theta / rad);
+ G4double wX = m_R[i] * sin(Theta) * cos(Phi);
+ G4double wY = m_R[i] * sin(Theta) * sin(Phi);
+ G4double wZ = m_R[i] * cos(Theta);
SuperX3w = G4ThreeVector(wX, wY, wZ);
// vector corresponding to the center of the module
SuperX3Center = SuperX3w;
// vector parallel to one axis of silicon plane
- G4double ii = cos(Theta / rad) * cos(Phi / rad);
- G4double jj = cos(Theta / rad) * sin(Phi / rad);
- G4double kk = -sin(Theta / rad);
+ G4double ii = cos(Theta) * cos(Phi);
+ G4double jj = cos(Theta) * sin(Phi);
+ G4double kk = -sin(Theta);
G4ThreeVector Y = G4ThreeVector(ii, jj, kk);
SuperX3w = SuperX3w.unit();
@@ -265,7 +265,7 @@ void SuperX3::ConstructDetector(G4LogicalVolume* world) {
SuperX3rot->rotate(m_beta_v[i], SuperX3v);
SuperX3rot->rotate(m_beta_w[i], SuperX3w);
// translation to place Telescope
- SuperX3pos = SuperX3w * Length * 0.5 + SuperX3Center;
+ SuperX3pos = SuperX3w * (SiliconThickness + AluStripThickness) * 0.5 + SuperX3Center;
}
VolumeMaker(i + 1, SuperX3pos, SuperX3rot, world);
diff --git a/Projects/SuperX3/SX3_20cm.detector b/Projects/SuperX3/SX3_20cm.detector
new file mode 100644
index 000000000..4f2220b26
--- /dev/null
+++ b/Projects/SuperX3/SX3_20cm.detector
@@ -0,0 +1,26 @@
+%%%%%%%
+Target
+ THICKNESS= 0.1 micrometer
+ ANGLE= 0 deg
+ RADIUS= 15 mm
+ MATERIAL= CD2
+ X= 0 mm
+ Y= 0 mm
+ Z= -20 cm
+ NbSlices= 10
+
+%%%%%%% Alias Phi
+% will create a new block for each value of Phi
+Alias Phi
+ Action= Copy
+ Value= 0 15 30 45 60 75 90 105 120 135 150 165 180 195 210 225 240 255 270 285 300 315 330 345
+
+%%%%%%% Detector 1 %%%%%%%
+SuperX3
+ THETA= 90 deg
+ PHI= @Phi deg
+ R= 198.4 mm
+ BETA= 0 0 0 deg
+ VIS= all
+
+
--
GitLab