deleted: Armory/ClassPC1An.h delete4d duplicate

modified:   Armory/ClassPW.h changes made for 4 wire shift to the anode and cathode starting positions
	modified:   Armory/PC_StepLadder_Correction.h
	modified:   anasen_fem/paraview_plotter.py corrected path for new pc
	modified:   run_17F.sh
	modified:   run_27Al.sh
	modified:   run_sx3.sh

Armory/ClassPC1An.h

@@ -1,246 +0,0 @@
-#ifndef ClassPC_h
-#define ClassPC_h
-
-#include <cstdio>
-#include <TMath.h>
-#include <TVector3.h>
-#include <TRandom.h>
-
-struct PCHit_1An{
-  std::pair<short, short>   nearestWire; // anode, cathode
-  std::pair<double, double> nearestDist; // anode, cathode
-
-  short   nextNearestWire; // cathode
-  double nextNearestDist; // cathode
-
-  void Clear(){
-    nearestWire.first = -1;
-    nearestWire.second = -1;
-    nearestDist.first = 999999999;
-    nearestDist.second = 999999999;
-    nextNearestWire= -1;
-    nextNearestDist = 999999999;
-  }
-};
-
-//!########################################################
-class PC{ // proportional wire
-public: 
-  PC(){ ClearHitInfo();};
-  ~PC(){};
-  
-  PCHit_1An GetHitInfo() const {return hitInfo;}
-  std::pair<short, short>   GetNearestID()          const {return hitInfo.nearestWire;}
-  std::pair<double, double> GetNearestDistance()    const {return hitInfo.nearestDist;}
-  short   Get2ndNearestID()       const {return hitInfo.nextNearestWire;}
-  double Get2ndNearestDistance() const {return hitInfo.nextNearestDist;}
-
-  TVector3 GetTrackPos() const {return trackPos;}
-  TVector3 GetTrackVec() const {return trackVec;}
-  double GetTrackTheta() const {return trackVec.Theta();}
-  double GetTrackPhi()   const {return trackVec.Phi();}
-  double GetZ0();
-
-  int GetNumWire() const {return nWire;}
-  double GetDeltaAngle() const {return dAngle;}
-  double GetAnodeLength() const {return anodeLength;}
-  double GetCathodeLength() const {return cathodeLength;}
-  TVector3 GetAnodeDn(short id) const {return An[id].first;}
-  TVector3 GetAnodeUp(short id) const {return An[id].second;}
-  TVector3 GetCathodeDn(short id) const {return Ca[id].first;}
-  TVector3 GetCathodeUp(short id) const {return Ca[id].second;}
-
-  TVector3 GetAnodneMid(short id)  const {return (An[id].first + An[id].second) * 0.5; }
-  double   GetAnodeTheta(short id) const {return (An[id].first - An[id].second).Theta();}
-  double   GetAnodePhi(short id)   const {return (An[id].first - An[id].second).Phi();}
-
-  TVector3 GetCathodneMid(short id)  const {return (Ca[id].first + Ca[id].second) * 0.5; }
-  double   GetCathodeTheta(short id) const {return (Ca[id].first - Ca[id].second).Theta();}
-  double   GetCathodePhi(short id)   const {return (Ca[id].first - Ca[id].second).Phi();}
-
-  void ClearHitInfo();
-  void ConstructGeo();
-  void FindWireID(TVector3 pos, TVector3 direction, bool verbose = false);
-  void CalTrack3(TVector3 sx3Pos, PCHit_1An hitInfo, double sigmaA = 0, double sigmaC = 0, bool verbose = false);
-
-  void Print(){
-    printf("     The nearest | Anode: %2d(%5.2f) Cathode: %2d(%5.2f)\n", hitInfo.nearestWire.first, 
-                                                                         hitInfo.nearestDist.first, 
-                                                                         hitInfo.nearestWire.second, 
-                                                                         hitInfo.nearestDist.second);
-
-    printf(" The 2nd nearest Cathode: %2d(%5.2f)\n", hitInfo.nextNearestWire, 
-                                                     hitInfo.nextNearestDist);
-  }
-
-private:
-
-  // PCHitInfo hitInfo;
-  PCHit_1An hitInfo;
-
-  TVector3 trackPos;
-  TVector3 trackVec;
-
-  const int nWire = 24;
-  const int wireShift = 3;
-  const float zLen = 380; //mm
-  const float radiusA = 37;
-  const float radiusC = 43;
-
-  double dAngle;
-  double anodeLength;
-  double cathodeLength;
-
-  std::vector<std::pair<TVector3,TVector3>> An; // the anode wire position vector in space 
-  std::vector<std::pair<TVector3,TVector3>> Ca; // the cathode wire position vector in space 
-
-  double Distance(TVector3 a1, TVector3 a2, TVector3 b1, TVector3 b2){
-    TVector3 na = a1 - a2;
-    TVector3 nb = b1 - b2;
-    TVector3 nd = (na.Cross(nb)).Unit();
-    return TMath::Abs(nd.Dot(a1-b2));
-  } 
-
-};
-
-inline void PC::ClearHitInfo(){
-  hitInfo.Clear();
-}
-
-inline void PC::ConstructGeo(){
-
-  An.clear();
-  Ca.clear();
-
-  std::pair<TVector3, TVector3> p1; // anode
-  std::pair<TVector3, TVector3> q1; // cathode
-
-  //anode and cathode start at pos-Y axis and count in right-Hand
-  //anode wire shift is right-hand.
-  //cathode wire shift is left-hand.
-
-  for(int i = 0; i < nWire; i++ ){
-    // Anode rotate right-hand
-    p1.first.SetXYZ( radiusA * TMath::Cos( TMath::TwoPi() / nWire * (i)  + TMath::PiOver2()),
-                    radiusA * TMath::Sin( TMath::TwoPi() / nWire * (i)  + TMath::PiOver2()),
-                    zLen/2);
-    p1.second.SetXYZ( radiusA * TMath::Cos( TMath::TwoPi() / nWire * (i + wireShift) + TMath::PiOver2()),
-                      radiusA * TMath::Sin( TMath::TwoPi() / nWire * (i + wireShift) + TMath::PiOver2()),
-                    -zLen/2);
-    An.push_back(p1);
-
-    // Cathod rotate left-hand
-    q1.first.SetXYZ( radiusC * TMath::Cos( TMath::TwoPi() / nWire * (i)  + TMath::PiOver2()),
-                    radiusC * TMath::Sin( TMath::TwoPi() / nWire * (i)  + TMath::PiOver2()),
-                    zLen/2);
-    q1.second.SetXYZ( radiusC * TMath::Cos( TMath::TwoPi() / nWire * (i - wireShift) + TMath::PiOver2()),
-                      radiusC * TMath::Sin( TMath::TwoPi() / nWire * (i - wireShift) + TMath::PiOver2()),
-                      -zLen/2);
-    Ca.push_back(q1);
-  }
-
-  dAngle = wireShift * TMath::TwoPi() / nWire;
-  anodeLength = TMath::Sqrt( zLen*zLen + TMath::Power(2* radiusA * TMath::Sin(dAngle/2),2) );
-  cathodeLength = TMath::Sqrt( zLen*zLen + TMath::Power(2* radiusC * TMath::Sin(dAngle/2),2) );
-}
-
-inline void PC::FindWireID(TVector3 pos, TVector3 direction, bool verbose ){
-
-  hitInfo.Clear();
-  double phi = direction.Phi();
-
-  for( int i = 0; i < nWire; i++){
-
-    double disA = 99999999;
-    double phiS = An[i].first.Phi()  - TMath::PiOver4();
-    double phiL = An[i].second.Phi() + TMath::PiOver4();
-    // printf("A%2d: %f %f | %f\n", i, phiS * TMath::RadToDeg(), phiL * TMath::RadToDeg(), phi * TMath::RadToDeg());
-    if( phi > 0 && phiS > phiL )  phiL = phiL + TMath::TwoPi();
-    if( phi < 0 && phiS > phiL )  phiS = phiS - TMath::TwoPi();
-
-    if( phiS < phi && phi < phiL) {
-      disA = Distance( pos, pos + direction, An[i].first, An[i].second);
-      if( disA < hitInfo.nearestDist.first ){
-        hitInfo.nearestDist.first = disA;
-        hitInfo.nearestWire.first = i;
-      }
-    }
-
-    double disC = 99999999;
-    phiS = Ca[i].second.Phi()- TMath::PiOver4();
-    phiL = Ca[i].first.Phi() + TMath::PiOver4();
-    // printf("C%2d: %f %f\n", i, phiS * TMath::RadToDeg(), phiL * TMath::RadToDeg());
-    if( phi > 0 && phiS > phiL ) phiL = phiL + TMath::TwoPi();
-    if( phi < 0 && phiS > phiL ) phiS = phiS - TMath::TwoPi();
-
-    if(phiS < phi && phi < phiL) {
-      disC = Distance( pos, pos + direction, Ca[i].first, Ca[i].second);
-      if( disC < hitInfo.nearestDist.second ){
-        hitInfo.nearestDist.second = disC;
-        hitInfo.nearestWire.second = i;
-      }
-    }
-
-    if(verbose) printf(" %2d | %8.2f, %8.2f\n", i, disA, disC);
-  }
-
-  short cathode1 = hitInfo.nearestWire.second;
-  short ccc1 = cathode1 - 1; if( ccc1 < 0 ) ccc1 += nWire;
-  short ccc2 = (cathode1 + 1) % nWire; 
-
-   double haha1 = Distance( pos, pos + direction, Ca[ccc1].first, Ca[ccc1].second);
-   double haha2 = Distance( pos, pos + direction, Ca[ccc2].first, Ca[ccc2].second);
-  if( haha1 < haha2){
-    hitInfo.nextNearestWire = ccc1;
-    hitInfo.nextNearestDist = haha1;
-  }else{
-    hitInfo.nextNearestWire = ccc2;
-    hitInfo.nextNearestDist= haha2;
-  }
-
-  if( verbose ) Print();
-}
-
-
-inline void PC::CalTrack3(TVector3 sx3Pos, PCHit_1An hitInfo, double sigmaA, double sigmaC, bool verbose){
-
-  trackPos = sx3Pos;
-
-  double p1 = TMath::Abs(hitInfo.nearestDist.first + gRandom->Gaus(0, sigmaA));
-  short anodeID1 = hitInfo.nearestWire.first;
-
-  double q1 = TMath::Abs(hitInfo.nearestDist.second + gRandom->Gaus(0, sigmaC));
-  double q2 = TMath::Abs(hitInfo.nextNearestDist+ gRandom->Gaus(0, sigmaC));
-  double fracC = q1 / (q1 + q2);
-  short cathodeID1 = hitInfo.nearestWire.second;
-  short cathodeID2 = hitInfo.nextNearestWire;
-  TVector3 shiftC1 = (Ca[cathodeID2].first - Ca[cathodeID1].first) * fracC;
-  TVector3 shiftC2 = (Ca[cathodeID2].second - Ca[cathodeID1].second) * fracC;
-
-  TVector3 a1 = An[anodeID1].first;
-
-  TVector3 c1 = Ca[cathodeID1].first + shiftC1;
-  TVector3 c2 = Ca[cathodeID1].second + shiftC2;
-
-  TVector3 n1 = (sx3Pos - a1).Unit();
-  TVector3 n2 = (c1 - c2).Cross((sx3Pos - c2)).Unit();
-
-  // if the handiness of anode and cathode revered, it should be n2 cross n1
-  trackVec = (n2.Cross(n1)).Unit();
-
-  if( verbose ) printf("Theta, Phi = %f, %f \n", trackVec.Theta() *TMath::RadToDeg(), trackVec.Phi()*TMath::RadToDeg()); 
-
-}
-
-inline double PC::GetZ0(){
-
-  double x = trackPos.X();
-  double y = trackPos.Y();
-  double rho = TMath::Sqrt(x*x + y*y);
-  double theta = trackVec.Theta();
-  
-  return trackPos.Z() - rho / TMath::Tan(theta);
-
-}
-
-#endif 

Armory/ClassPW.h

@@ -154,8 +154,8 @@ inline void PW::ConstructGeo()
   std::pair<TVector3, TVector3> q1; // cathode
 
   double k = TMath::TwoPi() / 24.; // 48 solder thru holes, wires in every other one
-  double offset_a1 = -6 * k - 3 * k;
+  double offset_a1 = -5 * k - 5 * k;
-  double offset_c1 = -4 * k - 2 * k - TMath::TwoPi() / 48; // correct for a half-turn
+  double offset_c1 = -5 * k - 2 * k - TMath::TwoPi() / 48; // correct for a half-turn
   // std::cerr << "Here!" << std::endl;
   // #include "../scratch/testing.h"
   double offset_a2 = offset_a1 + wireShift * k;

Armory/PC_StepLadder_Correction.h

@@ -32,7 +32,7 @@ double model_invert(double* y, double* p) {
  			break;
  		}
  	}
-    return result+80;
+    return result;
 }
 
 double model_a1c1(double* x, double* p) {

anasen_fem/paraview_plotter.py

@@ -1,5 +1,5 @@
-#!/home/vsitaraman/ParaView-6.1.0-RC1-MPI-Linux-Python3.12-x86_64/bin/pvbatch
+#!/home/vsitaraman/ParaView-6.1.0-MPI-Linux-Python3.12-x86_64/bin/pvbatch
-######## !/home/vs19g/ParaView-6.1.0-MPI-Linux-Python3.12-x86_64/bin/pvbatch
+########### !/home/vsitaraman/ParaView-6.1.0-RC1-MPI-Linux-Python3.12-x86_64/bin/pvbatch
 import numpy as np
 import sys
 from paraview.simple import *

run_17F.sh

@@ -2,7 +2,8 @@ rm results_run*.root
 export DATASET="17F"
 export flip180="0"
 export flipa=0
-export anode_offset=0
+export anode_offset=2
+export cathode_offset=0
 #root -q -l -b -x ../ANASEN_analysis/data/17F_Data/Source_005_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run05.root;
 #root -q -l -b -x ../ANASEN_analysis/data/17F_Data/Source_006_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run06.root;
 #root -q -l -b -x ../ANASEN_analysis/data/17F_Data/Source_007_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run07.root;
@@ -26,18 +27,19 @@ export anode_offset=0
 #root -q -l -b -x ../ANASEN_analysis/data/17F_Data/SourceRun_021_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run21.root;
 
 #17F reaction data
-export flip180="0"
+# export flip180="0"
-declare -i run=231 #49
+declare -i run=322 #49
-while [[ $run -lt 258 ]]; do #392
+while [[ $run -lt 399 ]]; do #392
     wrun=$(printf "%03d" $run)
     file_exists=$(test -f ../ANASEN_analysis/data/17F_Data/Run_"$wrun"_mapped.root)
     if [[ $file_exists -ne 0 ]]; then continue; fi
-    root -q -l -b -x ../ANASEN_analysis/data/17F_Data/Run_"$wrun"_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run$wrun.root;
+    root -q -l -b -x ../ANASEN_analysis/data/17F_Data/Run_"$wrun"_mapped.root -e 'tree->Process("MakeVertex.C+O","Analyzer_17F.root")'; mv Analyzer_17F.root 17F_output/results_run$wrun.root;
     run=run+1
 done
 
-rm output_17F.root
+rm output.root
-hadd -j 4 -k output_17F.root results_run2*.root
+hadd -k -j 4 output.root 17F_output/results_run*.root
+mv output.root output_17F.root
 
 unset souce_vertex
 unset DATASET

run_27Al.sh

@@ -2,7 +2,9 @@ rm results_run*.root
 export DATASET="27Al"
 export flip180="0"
 export flipa=0
+export flipc=0
 export anode_offset=0
+export cathode_offset=0
 #declare -i run=28
 #while [[ $run -lt 34 ]]; do #runs 1 to 84
 #    wrun=$(printf "%03d" $run)

run_sx3.sh

@@ -2,6 +2,8 @@
 # rm results_run*.root
 export flipa=0
 export anode_offset=0
+export cathode_offset=0
+
 export DATASET="27Al"
 if [[ 1 -eq 0 ]]; then
 #root -b -q -l -x ../ANASEN_analysis/data/27Al_Data/Run_009_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run09.root;
@@ -18,9 +20,8 @@ fi
 #exit
 #alpha+gas 27Al
 export DATASET="27Al"
-export flip180="0"
 #root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_009_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run09.root;
-if [[ 1 -eq 0 ]]; then
+if [[ 1 -eq 1 ]]; then
 #export timecut_low=230.0;
 export timecut_low=400.0;
 #export timecut_high=400.0;
@@ -29,7 +30,7 @@ export timecut_low=400.0;
 #export source_vertex=53.44; root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_010_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run10.root;
 #export source_vertex=53.44; root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_011_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run11.root;
 export source_vertex=53.44; root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_012_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run12.root;
-exit
+# exit
 #export source_vertex=53.44; root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_013_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run13.root;
 #exit
 fi
@@ -37,12 +38,10 @@ fi
 #protons+gas, 27Al
 #export flip180="1"
 #export flip180="0"
-if [[ 1 -eq 1 ]]; then
+if [[ 1 -eq 0 ]]; then
-export flipa=0
-export anode_offset=0
 export source_vertex=-200.0; #put the 'source' on the entrance window
 root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_018_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run18.root;
-exit
+# exit
 root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_015_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run15.root;
 root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_017_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run17.root;
 root -q -b -x ../ANASEN_analysis/data/27Al_Data/Run_019_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run19.root;
@@ -60,7 +59,6 @@ fi
 #root -l -x  results_run19.root results_run12.root -e "new TBrowser"
 #exit
 export DATASET="17F"
-export flip180="0"
 if [[ 1 -eq 0 ]]; then
 root -q -l -b -x ../ANASEN_analysis/data/17F_Data/Source_005_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run05.root;
 root -q -l -b -x ../ANASEN_analysis/data/17F_Data/Source_006_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run06.root;
@@ -86,8 +84,7 @@ export source_vertex=-24.96; root -q -l -b -x ../ANASEN_analysis/data/17F_Data/S
 export source_vertex=-73.96; root -q -l -b -x ../ANASEN_analysis/data/17F_Data/SourceRun_021_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run21.root;
 fi
 #17F reaction data
-#export flip180="0"
+if [[ 1 -eq 1 ]]; then
-if [[ 1 -eq 0 ]]; then
 export source_vertex=-57.28; root -q -l -b -x ../ANASEN_analysis/data/17F_Data/ProtonRun_035_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run35.root;
 #export source_vertex=-8.28; root -q -l -b -x ../ANASEN_analysis/data/17F_Data/ProtonRun_036_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root resulrs_run36.root;
 #export source_vertex=-27.88; root -q -l -b -x ../ANASEN_analysis/data/17F_Data/ProtonRun_037_mapped.root -e 'tree->Process("MakeVertex.C+O")'; mv Analyzer_SX3.root results_run37.root;