Commit 1131bddf authored by Vratislav Chudoba's avatar Vratislav Chudoba

Check of coordinates system finished. Script works at Ubunutu 18.04 with new EXPERTroot

parent 64f0f092
#if !defined(__CLING__)
#include "TSystem.h" #include "TSystem.h"
#include "TFile.h" #include "TFile.h"
#include "TTree.h" #include "TTree.h"
#include "TVector3.h" #include "TVector3.h"
//#include "../src/TNeEvent.h" #include "TH2F.h"
//#include "../../AculUtils/TELoss/TELoss.h" #include "../src/TNeEvent.h"
#include "../src/TNeDet16.h"
#include "TStopwatch.h"
#endif
using std::cout; using std::cout;
using std::endl; using std::endl;
...@@ -295,8 +301,9 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents ...@@ -295,8 +301,9 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents
//left telescope //left telescope
const Float_t z1mm = 230.; const Float_t z1mm = 230.;
const Float_t zThin = 230.-53.6; const Float_t zThin = 230.-53.6;
const Float_t xThinOffset = -3.; //TODO: check signs and values
const Float_t yThinOffset = -1.8; // const Float_t xThinOffset = -3.;
// const Float_t yThinOffset = -1.8;
//TODO: add sign taking into account direction of numbering //TODO: add sign taking into account direction of numbering
...@@ -304,11 +311,17 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents ...@@ -304,11 +311,17 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents
const Int_t kSQL_Y_strips = 16; const Int_t kSQL_Y_strips = 16;
const Int_t kSQL_20_strips = 16; const Int_t kSQL_20_strips = 16;
//1mm left size
const Float_t kSQL_X_size = 60.;
const Float_t kSQL_Y_size = 60.;
//thin detector //thin detector
const Double_t kSQ20_norm_thickness = 20.; const Double_t kSQ20_norm_thickness = 20.;
const Float_t kSQ20_X_size = 50.;
const Float_t kSQ20_Y_size = 50.;
const Double_t thinXoffset = 1.; const Double_t kThinXoffset = -1.;
const Double_t thinYoffset = -1.8; const Double_t kThinYoffset = -1.8;
//left CsI detectors //left CsI detectors
...@@ -672,6 +685,9 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents ...@@ -672,6 +685,9 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents
// correction for thickness inhomogeneity // correction for thickness inhomogeneity
// when not taking into account hitting // when not taking into account hitting
// particle tracking // particle tracking
//
// note: correction for thickness is not affected
// by used system of coordinates
/////////////////////////////////////////////// ///////////////////////////////////////////////
Double_t currThickness; Double_t currThickness;
...@@ -679,13 +695,14 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents ...@@ -679,13 +695,14 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents
if (SQLXmult == 1 && SQLYmult == 1) { if (SQLXmult == 1 && SQLYmult == 1) {
for (Int_t yi = 0; yi < kSQL_Y_strips; yi++) { for (Int_t yi = 0; yi < kSQL_Y_strips; yi++) {
// cout << "y bin: " << yi+1 << "\t\t";
//if energy in 1mm Y strip is reasonable
if (SQLYE[yi]>0) { if (SQLYE[yi]>0) {
for (Int_t xi = 0; xi < TMath::Abs(kSQL_X_strips); xi++) { for (Int_t xi = 0; xi < TMath::Abs(kSQL_X_strips); xi++) {
//if energy in 1mm X strip is reasonable
if (SQLXE[xi] > 0) { if (SQLXE[xi] > 0) {
//TODO: check the mapping here at first
currThickness = hThickness->GetBinContent(xi+1, yi+1); currThickness = hThickness->GetBinContent(xi+1, yi+1);
for (Int_t x20 = 0; x20 < kSQL_20_strips; x20++) { for (Int_t x20 = 0; x20 < kSQL_20_strips; x20++) {
...@@ -695,12 +712,12 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents ...@@ -695,12 +712,12 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents
SQ20EcorrMult++; SQ20EcorrMult++;
} }
}//for*/ }//for thin detector X strips
} }//if 1 mm X
} }//for thick detector X strips
} }//if 1 mm Y
} }//for thick detector Y strips
}//if for correction according energy multiplicity }//if for correction according energy multiplicity
...@@ -724,8 +741,8 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents ...@@ -724,8 +741,8 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents
if (SQLXtimeMult==1) { if (SQLXtimeMult==1) {
for (Int_t i = 0; i<32; i++) { for (Int_t i = 0; i<32; i++) {
if (SQLXEtimeFiltered[i]>0) { if (SQLXEtimeFiltered[i]>0) {
//todo: check the system of coordinates // x1mm = -30. + (i+1./2.)*60./32.;
x1mm = -30. + (i+1./2.)*60./32.; x1mm = (i+0.5)*kSQL_X_size/kSQL_X_strips + kSQL_X_size/2.; //ok
} }
} }
} }
...@@ -734,14 +751,15 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents ...@@ -734,14 +751,15 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents
if (SQLYtimeMult==1) { if (SQLYtimeMult==1) {
for (Int_t i = 0; i<16; i++) { for (Int_t i = 0; i<16; i++) {
if (SQLYEtimeFiltered[i]>0) { if (SQLYEtimeFiltered[i]>0) {
//todo: check the system of coordinates // y1mm = -30. + (i+1./2.)*60./16.;
y1mm = -30. + (i+1./2.)*60./16.; y1mm = -kSQL_Y_size/2. + (i+0.5)*kSQL_Y_size/kSQL_Y_strips; //ok
} }
} }
} }
//for MWPC wire multiplicity == 1 //for MWPC wire multiplicity == 1
//vHit1mm conversion to plane of thin detector
if (SQLXtimeMult==1 && SQLYtimeMult==1 && xt>-50. && yt>-50.) { if (SQLXtimeMult==1 && SQLYtimeMult==1 && xt>-50. && yt>-50.) {
vHit1mm.SetXYZ(x1mm-xt, y1mm-yt, z1mm); vHit1mm.SetXYZ(x1mm-xt, y1mm-yt, z1mm);
vHit1mm *= zThin/z1mm; vHit1mm *= zThin/z1mm;
...@@ -750,16 +768,22 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents ...@@ -750,16 +768,22 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents
if (SQLXtimeMult==1 && SQLYtimeMult==1 if (SQLXtimeMult==1 && SQLYtimeMult==1
&& vHit1mm.X()+xt>-30. && vHit1mm.X()+xt<30. //check if the particle passed through sensitive area of thin detector
&& vHit1mm.Y()+yt>-30. && vHit1mm.Y()+yt<30.) { && vHit1mm.X()+xt < (kSQ20_X_size/2. + kThinXoffset)
&& vHit1mm.X()+xt > (-kSQ20_X_size/2. + kThinXoffset)
&& vHit1mm.Y()+yt < (kSQ20_Y_size/2. + kThinYoffset)
&& vHit1mm.Y()+yt > (-kSQ20_Y_size/2. + kThinYoffset)
) {
mapYbin = (vHit1mm.Y()+yt+thinYoffset+30.-0.5*60./16.)*16./60.;
mapXbin = (vHit1mm.X()+xt+thinXoffset+30.-0.5*60./32.)*32./60.;
//coordinate in the thin detector plane //coordinate in the thin detector plane
xThin = vHit1mm.X()+xt; xThin = vHit1mm.X()+xt;
yThin = vHit1mm.Y()+yt; yThin = vHit1mm.Y()+yt;
mapYbin = yThin*kSQL_Y_strips/kSQL_Y_size + kSQL_Y_strips/2.; //ok
mapXbin = xThin*kSQL_X_strips/kSQL_X_size - kSQL_X_strips/2.; //ok
//calculation of corrected energy //calculation of corrected energy
//TODO: check the mapping here at second //TODO: check the mapping here at second
currThickness = hThickness->GetBinContent(mapXbin+1, mapYbin+1); currThickness = hThickness->GetBinContent(mapXbin+1, mapYbin+1);
...@@ -767,6 +791,7 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents ...@@ -767,6 +791,7 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents
for (Int_t x20 = 0; x20 < kSQL_20_strips; x20++) { for (Int_t x20 = 0; x20 < kSQL_20_strips; x20++) {
if (SQ20E[x20] > 0/*SQ20timeMult==1*/ && currThickness<30.) { if (SQ20E[x20] > 0/*SQ20timeMult==1*/ && currThickness<30.) {
SQ20EcorrHit[x20] = (kSQ20_norm_thickness/currThickness)*SQ20E[x20]; SQ20EcorrHit[x20] = (kSQ20_norm_thickness/currThickness)*SQ20E[x20];
//when multiplicity larger than 1, we obtain incorrect information
SQ20EcorrHitMult++; SQ20EcorrHitMult++;
// cout << mapXbin << "\t" << mapYbin << "\t" << x20 // cout << mapXbin << "\t" << mapYbin << "\t" << x20
// << "\t" << currThickness // << "\t" << currThickness
...@@ -789,16 +814,26 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents ...@@ -789,16 +814,26 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents
if (SQLXtimeMult==1 && SQLYtimeMult==1 if (SQLXtimeMult==1 && SQLYtimeMult==1
&& vHit1mm.X()+xtc>-30. && vHit1mm.X()+xtc<30. //check if the particle passed through sensitive area of thin detector
&& vHit1mm.Y()+ytc>-30. && vHit1mm.Y()+ytc<30.) { && vHit1mm.X()+xtc < (kSQ20_X_size/2. + kThinXoffset)
&& vHit1mm.X()+xtc > (-kSQ20_X_size/2. + kThinXoffset)
&& vHit1mm.Y()+ytc < (kSQ20_Y_size/2. + kThinYoffset)
&& vHit1mm.Y()+ytc > (-kSQ20_Y_size/2. + kThinYoffset)
// && vHit1mm.X()+xtc>-30. && vHit1mm.X()+xtc<30.
// && vHit1mm.Y()+ytc>-30. && vHit1mm.Y()+ytc<30.
) {
mapYbin = (vHit1mm.Y()+ytc+thinYoffset+30.-0.5*60./16.)*16./60.;
mapXbin = (vHit1mm.X()+xtc+thinXoffset+30.-0.5*60./32.)*32./60.;
//coordinate in the thin detector plane
xThin = vHit1mm.X()+xtc; xThin = vHit1mm.X()+xtc;
yThin = vHit1mm.Y()+ytc; yThin = vHit1mm.Y()+ytc;
mapYbin = yThin*kSQL_Y_strips/kSQL_Y_size + kSQL_Y_strips/2.; //ok
mapXbin = xThin*kSQL_X_strips/kSQL_X_size - kSQL_X_strips/2.; //ok
// mapYbin = (vHit1mm.Y()+ytc+kThinYoffset+30.-0.5*60./16.)*16./60.;
// mapXbin = (vHit1mm.X()+xtc+kThinXoffset+30.-0.5*60./32.)*32./60.;
//coordinate in the thin detector plane
//calculation of corrected energy //calculation of corrected energy
currThickness = hThickness->GetBinContent(mapXbin+1, mapYbin+1); currThickness = hThickness->GetBinContent(mapXbin+1, mapYbin+1);
...@@ -814,7 +849,8 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents ...@@ -814,7 +849,8 @@ void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents
}//for*/ }//for*/
if (SQ20EcorrHitCMult==1) { if (SQ20EcorrHitCMult==1) {
expectedThinStrip = (xThin + 25. - thinXoffset)*(16./50.)-0.5; //TODO: check this statement
expectedThinStrip = (xThin + 25. - kThinXoffset)*(16./50.)-0.5;
} }
}//if inside the map }//if inside the map
...@@ -849,6 +885,8 @@ void fillChain(const TString beam = "he", const Int_t from = 0, const Int_t to = ...@@ -849,6 +885,8 @@ void fillChain(const TString beam = "he", const Int_t from = 0, const Int_t to =
TStopwatch stopwatch; TStopwatch stopwatch;
stopwatch.Start(); stopwatch.Start();
// cout << "aksjdlkajbdlkjasd" << endl;
for (Int_t i = from; i <= to; i++) { for (Int_t i = from; i <= to; i++) {
fillTree(beam, i, noevents); fillTree(beam, i, noevents);
} }
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment