diff --git a/.vscode/settings.json b/.vscode/settings.json index e1fb437..8ae7276 100644 --- a/.vscode/settings.json +++ b/.vscode/settings.json @@ -15,6 +15,7 @@ "ApplyMapping.C": "cpp", "script.C": "cpp", "Analyzer.C": "cpp", - "PreAnalysis.C": "cpp" + "PreAnalysis.C": "cpp", + "ANASEN_model.C": "cpp" } } \ No newline at end of file diff --git a/ANASEN_model.C b/ANASEN_model.C index b873f73..ee4a3cf 100644 --- a/ANASEN_model.C +++ b/ANASEN_model.C @@ -61,7 +61,7 @@ void ANASEN_model(int anodeID1 = -1, int anodeID2 = -1, int cathodeID1 = -1, int pcA->SetLineColor(4); for( int i = 0; i < nWire; i++){ - if( i < anodeID1 || i > anodeID2) continue; + if( anodeID2 >= 0 && (i < anodeID1 || i > anodeID2) ) continue; worldBox->AddNode(pcA, i+1, new TGeoCombiTrans( radiusAnew * TMath::Cos( TMath::TwoPi() / nWire *i + dAngle / 2), radiusAnew * TMath::Sin( TMath::TwoPi() / nWire *i + dAngle / 2), 0, @@ -76,7 +76,7 @@ void ANASEN_model(int anodeID1 = -1, int anodeID2 = -1, int cathodeID1 = -1, int TGeoVolume *pcC = geom->MakeTube("tub2", Al, 0, 0.01, wireCLength/2); pcC->SetLineColor(6); for( int i = 0; i < nWire; i++){ - if( i < cathodeID1 || i > cathodeID2) continue; + if( cathodeID2 >= 0 && (i < cathodeID1 || i > cathodeID2) ) continue; worldBox->AddNode(pcC, i+1, new TGeoCombiTrans( radiusCnew * TMath::Cos( TMath::TwoPi() / nWire *i - dAngle/2), radiusCnew * TMath::Sin( TMath::TwoPi() / nWire *i - dAngle/2), 0, @@ -92,15 +92,26 @@ void ANASEN_model(int anodeID1 = -1, int anodeID2 = -1, int cathodeID1 = -1, int TGeoVolume * sx3 = geom->MakeBox("box", Al, 0.1, sx3Width/2, sx3Length/2); sx3->SetLineColor(kGreen+3); for( int i = 0; i < nSX3; i++){ - worldBox->AddNode(sx3, 2*i+1., new TGeoCombiTrans( sx3Radius * TMath::Cos( TMath::TwoPi() / nSX3 *i), - sx3Radius * TMath::Sin( TMath::TwoPi() / nSX3 *i), + worldBox->AddNode(sx3, 2*i+1., new TGeoCombiTrans( sx3Radius * TMath::Cos( TMath::TwoPi() / nSX3 * (i + 0.5)), + sx3Radius * TMath::Sin( TMath::TwoPi() / nSX3 * (i + 0.5)), sx3Length/2+sx3Gap, - new TGeoRotation("rot1", 360/nSX3 * (i), 0., 0.))); + new TGeoRotation("rot1", 360/nSX3 * (i + 0.5), 0., 0.))); - worldBox->AddNode(sx3, 2*i+2., new TGeoCombiTrans( sx3Radius* TMath::Cos( TMath::TwoPi() / nSX3 *i), - sx3Radius * TMath::Sin( TMath::TwoPi() / nSX3 *i), - -sx3Length/2-sx3Gap, - new TGeoRotation("rot1", 360/nSX3 * (i), 0., 0.))); + worldBox->AddNode(sx3, 2*i+2., new TGeoCombiTrans( sx3Radius * TMath::Cos( TMath::TwoPi() / nSX3 * (i + 0.5)), + sx3Radius * TMath::Sin( TMath::TwoPi() / nSX3 * (i + 0.5)), + -sx3Length/2-sx3Gap, + new TGeoRotation("rot1", 360/nSX3 * (i + 0.5), 0., 0.))); + } + + const int qqqR1 = 10; + const int qqqR2 = 50; + TGeoVolume *qqq = geom->MakeTubs("qqq", Al, qqqR1, qqqR2, 0.5, 5, 85); + qqq->SetLineColor(7); + for( int i = 0; i < 4; i++){ + worldBox->AddNode(qqq, i+1, new TGeoCombiTrans( 0, + 0, + 100, + new TGeoRotation("rot1", 360/4 * (i), 0., 0.))); } geom->CloseGeometry(); diff --git a/Analyzer.C b/Analyzer.C index 3b47c05..be89e9c 100644 --- a/Analyzer.C +++ b/Analyzer.C @@ -126,6 +126,8 @@ Bool_t Analyzer::Process(Long64_t entry){ //########################################################### Track constrcution + //############################## DO THE KINEMATICS + return kTRUE; } diff --git a/ClassAnasen.h b/ClassAnasen.h index a94a02e..64873ae 100644 --- a/ClassAnasen.h +++ b/ClassAnasen.h @@ -5,10 +5,19 @@ #include #include -class Anasen{ +#include "TGeoManager.h" +#include "TGeoVolume.h" +#include "TGeoBBox.h" +#include "TCanvas.h" +#include "TPolyMarker3D.h" +#include "TPolyLine3D.h" + +class ANASEN{ public: - Anasen(); - ~Anasen() {} + ANASEN(); + ~ANASEN(); + + TVector3 CalSX3Pos(unsigned short ID, unsigned short chUp, unsigned short chDown, unsigned short chBack, float eUp = 0, float eDown = 0 ); void CalTrack(TVector3 sx3Pos, int anodeID, int cathodeID); @@ -18,22 +27,25 @@ public: double GetTrackTheta() const {return trackVec.Theta();} double GetTrackPhi() const {return trackVec.Phi();} - // void DrawTrack(); + void DrawAnasen(int anodeID1 = -1, int anodeID2 = -1, int cathodeID1 = -1, int cathodeID2 = -1, bool DrawQQQ = false ); private: const int nWire = 24; - const int wireShift = 3; // how twisted is the wire. - const double AnodeRadius = 38 ; // mm - const double CathodeRadius = 38 ; // mm - const double PCLength = 100; //mm + const int wireShift = 3; + const int zLen = 380; //mm + const int radiusA = 37; + const int radiusC = 43; std::vector P1; // the anode wire position vector in space std::vector P2; // the anode wire position vector in space std::vector Q1; // the cathode wire position vector in space std::vector Q2; // the cathode wire position vector in space - void CalWireDirection(); + std::vector> S1; // coners of the SX3 0-11, z = mid point + std::vector> S2; // coners of the SX3 12-23, z = mid point + + void CalGeometry(); TVector3 trackPos; TVector3 trackVec; @@ -44,15 +56,40 @@ private: TVector3 trackVecErrorA; // error vector prependicular to the Anode-Pos plan TVector3 trackVecErrorC; // error vector prependicular to the Cathode-Pos plan + const int nSX3 = 12; + const int sx3Radius = 88; + const int sx3Width = 40; + const int sx3Length = 75; + const int sx3Gap = 46; + + const int qqqR1 = 10; + const int qqqR2 = 50; + const int qqqZPos = sx3Gap/2 + sx3Length + 30; + + TGeoManager *geom; + TGeoVolume *worldBox; + + void Construct3DModel(int anodeID1 = -1, int anodeID2 = -1, int cathodeID1 = -1, int cathodeID2 = -1, bool DrawQQQ = true); + }; -inline Anasen::Anasen(){ +//============================================== +inline ANASEN::ANASEN(){ - CalWireDirection(); + CalGeometry(); + + geom = nullptr; + worldBox = nullptr; } -inline void Anasen::CalWireDirection(){ +inline ANASEN::~ANASEN(){ + + delete geom; + +} + +inline void ANASEN::CalGeometry(){ TVector3 p1; // anode TVector3 p2; @@ -61,31 +98,155 @@ inline void Anasen::CalWireDirection(){ for(int i = 0; i < nWire; i++ ){ // Anode rotate right-hand - - p1.SetXYZ( AnodeRadius * TMath::Cos( TMath::TwoPi() * i / nWire ), - AnodeRadius * TMath::Sin( TMath::TwoPi() * i / nWire ), - -PCLength/2); - p2.SetXYZ( AnodeRadius * TMath::Cos( TMath::TwoPi() * (i + wireShift) / nWire ), - AnodeRadius * TMath::Sin( TMath::TwoPi() * (i + wireShift) / nWire ), - PCLength/2); + p1.SetXYZ( radiusA * TMath::Cos( TMath::TwoPi() * i / nWire ), + radiusA * TMath::Sin( TMath::TwoPi() * i / nWire ), + zLen/2); + p2.SetXYZ( radiusA * TMath::Cos( TMath::TwoPi() * (i + wireShift) / nWire ), + radiusA * TMath::Sin( TMath::TwoPi() * (i + wireShift) / nWire ), + -zLen/2); P1.push_back(p1); P2.push_back(p2); // Cathod rotate left-hand - - p1.SetXYZ( CathodeRadius * TMath::Cos( TMath::TwoPi() * i / nWire ), - CathodeRadius * TMath::Sin( TMath::TwoPi() * i / nWire ), - -PCLength/2); - p2.SetXYZ( CathodeRadius * TMath::Cos( TMath::TwoPi() * (i - wireShift) / nWire ), - CathodeRadius * TMath::Sin( TMath::TwoPi() * (i - wireShift) / nWire ), - PCLength/2); + p1.SetXYZ( radiusC * TMath::Cos( TMath::TwoPi() * i / nWire ), + radiusC * TMath::Sin( TMath::TwoPi() * i / nWire ), + zLen/2); + p2.SetXYZ( radiusC * TMath::Cos( TMath::TwoPi() * (i - wireShift) / nWire ), + radiusC * TMath::Sin( TMath::TwoPi() * (i - wireShift) / nWire ), + -zLen/2); Q1.push_back(p1); Q2.push_back(p2); } + TVector3 sa, sb; + for(int i = 0; i < nSX3; i++){ + sa.SetXYZ( sx3Radius, sx3Width/2, sx3Gap/2 + sx3Length/2 ); + sb.SetXYZ( sx3Radius, -sx3Width/2, sx3Gap/2 + sx3Length/2 ); + + sa.RotateZ( TMath::TwoPi() / nSX3 * (i + 0.5) ); + sb.RotateZ( TMath::TwoPi() / nSX3 * (i + 0.5) ); + S1.push_back(std::pair(sa,sb)); + + sa.SetXYZ( sx3Radius, sx3Width/2, -sx3Gap/2 - sx3Length/2 ); + sb.SetXYZ( sx3Radius, -sx3Width/2, -sx3Gap/2 - sx3Length/2 ); + + sa.RotateZ( TMath::TwoPi() / nSX3 * (i + 0.5) ); + sb.RotateZ( TMath::TwoPi() / nSX3 * (i + 0.5) ); + S2.push_back(std::pair(sa,sb)); + } + } -inline void Anasen::CalTrack(TVector3 sx3Pos, int anodeID, int cathodeID){ +inline void ANASEN::Construct3DModel(int anodeID1, int anodeID2, int cathodeID1, int cathodeID2, bool DrawQQQ ){ + + if( geom ) delete geom; + + // Create ROOT manager and master volume + geom = new TGeoManager("Detector", "ANASEN"); + + //--- define some materials + TGeoMaterial *matVacuum = new TGeoMaterial("Vacuum", 0,0,0); + TGeoMaterial *matAl = new TGeoMaterial("Al", 26.98,13,2.7); + //--- define some media + TGeoMedium *Vacuum = new TGeoMedium("Vacuum",1, matVacuum); + TGeoMedium *Al = new TGeoMedium("Root Material",2, matAl); + + //--- make the top container volume + Double_t worldx = 200.; //mm + Double_t worldy = 200.; //mm + Double_t worldz = 200.; //mm + worldBox = geom->MakeBox("ROOT", Vacuum, worldx, worldy, worldz); + geom->SetTopVolume(worldBox); + + //--- making axis + TGeoVolume *axisX = geom->MakeTube("axisX", Al, 0, 0.1, 5.); + axisX->SetLineColor(1); + worldBox->AddNode(axisX, 1, new TGeoCombiTrans(5, 0, 0., new TGeoRotation("rotA", 90., 90., 0.))); + + TGeoVolume *axisY = geom->MakeTube("axisY", Al, 0, 0.1, 5.); + axisY->SetLineColor(1); + worldBox->AddNode(axisY, 1, new TGeoCombiTrans(0, 5, 0., new TGeoRotation("rotB", 0., 90., 0.))); + + TGeoVolume *axisZ = geom->MakeTube("axisZ", Al, 0, 0.1, 5.); + axisZ->SetLineColor(1); + worldBox->AddNode(axisZ, 1, new TGeoTranslation(0, 0, 5)); + + //.......... convert to wire center dimensions + double dAngle = wireShift * TMath::TwoPi() / nWire; + double radiusAnew = radiusA * TMath::Cos( dAngle / 2.); + double wireALength = TMath::Sqrt( zLen*zLen + TMath::Power(2* radiusA * TMath::Sin(dAngle/2),2) ); + double wireATheta = TMath::ATan2( 2* radiusA * TMath::Sin( dAngle / 2.), zLen); + + // printf(" dAngle : %f\n", dAngle); + // printf(" newRadius : %f\n", radiusAnew); + // printf("wireLength : %f\n", wireALength); + // printf("wire Theta : %f\n", wireATheta); + + TGeoVolume *pcA = geom->MakeTube("tub1", Al, 0, 0.01, wireALength/2); + pcA->SetLineColor(4); + + for( int i = 0; i < nWire; i++){ + if( anodeID2 >= 0 && (i < anodeID1 || i > anodeID2) ) continue; + worldBox->AddNode(pcA, i+1, new TGeoCombiTrans( radiusAnew * TMath::Cos( TMath::TwoPi() / nWire *i + dAngle / 2), + radiusAnew * TMath::Sin( TMath::TwoPi() / nWire *i + dAngle / 2), + 0, + new TGeoRotation("rot1", 360/ nWire * (i + wireShift/2.), wireATheta * 180/ TMath::Pi(), 0.))); + } + + double radiusCnew = radiusC * TMath::Cos( dAngle / 2.); + double wireCLength = TMath::Sqrt( zLen*zLen + TMath::Power(2* radiusC * TMath::Sin(dAngle/2),2) ); + double wireCTheta = TMath::ATan2( 2* radiusC * TMath::Sin( dAngle / 2.), zLen); + + TGeoVolume *pcC = geom->MakeTube("tub2", Al, 0, 0.01, wireCLength/2); + pcC->SetLineColor(6); + for( int i = 0; i < nWire; i++){ + if( cathodeID2 >= 0 && (i < cathodeID1 || i > cathodeID2) ) continue; + worldBox->AddNode(pcC, i+1, new TGeoCombiTrans( radiusCnew * TMath::Cos( TMath::TwoPi() / nWire *i - dAngle/2), + radiusCnew * TMath::Sin( TMath::TwoPi() / nWire *i - dAngle/2), + 0, + new TGeoRotation("rot1", 360/ nWire * (i - wireShift/2.), -wireCTheta * 180/ TMath::Pi(), 0.))); + } + + TGeoVolume * sx3 = geom->MakeBox("box", Al, 0.1, sx3Width/2, sx3Length/2); + sx3->SetLineColor(kGreen+3); + for( int i = 0; i < nSX3; i++){ + worldBox->AddNode(sx3, 2*i+1., new TGeoCombiTrans( sx3Radius * TMath::Cos( TMath::TwoPi() / nSX3 * (i + 0.5)), + sx3Radius * TMath::Sin( TMath::TwoPi() / nSX3 * (i + 0.5)), + sx3Length/2+sx3Gap/2, + new TGeoRotation("rot1", 360/nSX3 * (i + 0.5), 0., 0.))); + + worldBox->AddNode(sx3, 2*i+2., new TGeoCombiTrans( sx3Radius * TMath::Cos( TMath::TwoPi() / nSX3 * (i + 0.5)), + sx3Radius * TMath::Sin( TMath::TwoPi() / nSX3 * (i + 0.5)), + -sx3Length/2-sx3Gap/2, + new TGeoRotation("rot1", 360/nSX3 * (i + 0.5), 0., 0.))); + } + + if( DrawQQQ ){ + TGeoVolume *qqq = geom->MakeTubs("qqq", Al, qqqR1, qqqR2, 0.5, 5, 85); + qqq->SetLineColor(7); + for( int i = 0; i < 4; i++){ + worldBox->AddNode(qqq, i+1, new TGeoCombiTrans( 0, + 0, + qqqZPos, + new TGeoRotation("rot1", 360/4 * (i), 0., 0.))); + } + } + +} + + +inline void ANASEN::DrawAnasen(int anodeID1, int anodeID2, int cathodeID1, int cathodeID2, bool DrawQQQ ){ + + + Construct3DModel(anodeID1, anodeID2, cathodeID1, cathodeID2, DrawQQQ); + + geom->CloseGeometry(); + geom->SetVisLevel(4); + worldBox->Draw("ogle"); + +} + +inline void ANASEN::CalTrack(TVector3 sx3Pos, int anodeID, int cathodeID){ trackPos = sx3Pos; @@ -95,7 +256,42 @@ inline void Anasen::CalTrack(TVector3 sx3Pos, int anodeID, int cathodeID){ // if the handiness of anode and cathode revered, it should be n2 cross n1 trackVec = (n1.Cross(n2)).Unit(); +} +inline TVector3 ANASEN::CalSX3Pos(unsigned short ID, unsigned short chUp, unsigned short chDown, unsigned short chBack, float eUp, float eDown){ + + TVector3 haha; + + if( (chUp - chDown) != 1 || (chDown % 2) != 0) return haha; + + int reducedID = ID % nSX3; + + TVector3 sa, sb; + + if( ID < nSX3 ){ //down + + sa = S1[reducedID].first; + sb = S1[reducedID].second; + + }else{ + + sa = S2[reducedID].first; + sb = S2[reducedID].second; + + } + + haha.SetX( (sb.X() - sa.X()) * chUp/8 + sa.X()); + haha.SetY( (sb.Y() - sa.Y()) * chUp/8 + sa.Y()); + + if( eUp == 0 || eDown == 0 ){ + haha.SetZ( sa.Z() + (2*(chBack - 7)-1) * sx3Length / 8 ); + }else{ + double frac = (eUp - eDown)/(eUp + eDown); // from +1 (downstream) to -1 (upstream) + double zPos = sa.Z() + sx3Length * frac/2; + haha.SetZ( zPos ); + } + + return haha; } diff --git a/ClassDet.h b/ClassDet.h index 4b57ed3..7cb211a 100644 --- a/ClassDet.h +++ b/ClassDet.h @@ -15,7 +15,7 @@ public: unsigned short e[MAXMULTI]; unsigned long long t[MAXMULTI]; - unsigned short index[MAXMULTI]; // id * 12 + ch; + unsigned short index[MAXMULTI]; // id * nCh + ch; bool used[MAXMULTI]; void Clear(){ @@ -56,4 +56,4 @@ private: }; -#endif \ No newline at end of file +#endif