Subversion Repositories f9daq

Rev

Rev 72 | Go to most recent revision | Blame | Compare with Previous | Last modification | View Log | RSS feed

  1. #ifndef CVODNIK_H
  2. #define CVODNIK_H
  3.  
  4. #include "TROOT.h"
  5. #include "TVector3.h"
  6. #include "TMath.h"
  7. #include "TPolyLine3D.h"
  8. #include "TRandom.h"
  9. #include "TDatime.h"
  10. #include "TH1F.h"
  11. #include "TH2F.h"
  12. #include "TColor.h"
  13.  
  14.  
  15.  
  16. #define MARGIN 1.0e-6
  17. #define DEGREE 180.0/3.14159265358979312
  18. #define MAX_REFLECTIONS 19
  19.  
  20.  
  21. using namespace TMath;
  22.  
  23. const double c_reflectivity = 0.96;
  24. const int REFRACTION = 1;
  25. const int REFLECTION = 2;
  26. const int TOTAL = 3;
  27. const int ABSORBED = -2;
  28. const TVector3 CENTER(-2.0, 0.0, 0.0);
  29. const int nch = 50; // number of "channels" used as number of bins in histograms
  30.  
  31. enum Fate {missed=0, hitExit=1, refracted=-1, enter=2, noreflection=-2, backreflected=-3, rays=3, absorbed=4};
  32.  
  33. int dbg = 0;
  34.  
  35. //=================================================================================
  36. // zarek (premica)
  37. class CRay
  38. {
  39. public:
  40.   CRay() :
  41.     r(TVector3(0,0,0)),
  42.     k(TVector3(1,0,0)),
  43.     //p(TVector3(0,1,0)),
  44.     p(TVector3(0,0,1)),
  45.     color(kBlack)
  46.   {};
  47.   CRay(TVector3 r0, TVector3 k0) :
  48.     r(r0),
  49.     k(k0.Unit()),
  50.     //p(TVector3(0,1,0)),
  51.     p(TVector3(0,0,1)),
  52.     color(kBlack)
  53.   {};
  54.   CRay(double x0, double y0, double z0, double l0, double m0, double n0) :
  55.     r(TVector3(x0,y0,z0)),
  56.     k(TVector3(l0,m0,n0).Unit()),
  57.     //p(TVector3(0,1,0)),
  58.     p(TVector3(0,0,1)),
  59.     color(kBlack)
  60.   {};
  61.  
  62.   void Set(TVector3 r0, TVector3 k0);
  63.   //void Set(double x0, double y0, double z0, double l0, double m0, double n0);
  64.   void SetColor(int c){color = c;};
  65.   void setPolarization(TVector3 p0) {
  66.     p = p0.Unit();
  67.     if (p.Dot(k)>1e-3) printf("*** ERROR in CRay: E has component || with k\n");
  68.   };
  69.  
  70.   //inline CRay & operator = (const CRay &);
  71.  
  72.   TVector3 GetR() const {return r;};
  73.   TVector3 GetK() const {return k;};
  74.   TVector3 GetP() const {return p;};
  75.  
  76.   void Print();
  77.   void Draw();
  78.   void Draw(double x_from, double x_to);
  79.   void DrawS(double x_from, double t);
  80.  
  81.   //r = position, k = unit wave vector, p = polarization
  82. private:
  83.   TVector3 r;
  84.   TVector3 k;
  85.   TVector3 p;
  86.   int color;
  87. };
  88.  
  89.  
  90. //=================================================================================
  91. // ravnina, definirana in omejena s 4 tockami na njej (stirikotni poligon)
  92. // zanasam se na taksno definicijo: (normala kaze noter)
  93. // r1----->r2
  94. // ^        |
  95. // |        |
  96. // |        |
  97. // |
  98. // r0<-----r3
  99. class CPlane4
  100. {
  101.  
  102. public:
  103.   CPlane4();
  104.   CPlane4(TVector3 r1, TVector3 r2, TVector3 r3, TVector3 r4);
  105.   CPlane4(TVector3 *vr);
  106.   void Set(TVector3 r1, TVector3 r2, TVector3 r3, TVector3 r4);
  107.   void Set(TVector3 *vr) {Set(vr[0], vr[1], vr[2], vr[3]);};
  108.   void FlipN(){n = -n;};
  109.  
  110.   int GetIntersection(TVector3 *vec, CRay ray);
  111.   int IsInTri(TVector3 vec, TVector3 e1, TVector3 e2, double);
  112.   int IsVectorIn(TVector3 vec);
  113.   int TestIntersection(CRay in); // ray in & ??? plane
  114.   int TestIntersection(TVector3 *vec, CRay in); // plane defined with normal vec and ray in
  115.  
  116.   TVector3 GetN() {return n;};
  117.  
  118.   void Print();
  119.   void Draw(int color = 1, int width = 1);
  120.  
  121. private:
  122.   TVector3 r[4], n;
  123.   double A, B, C, D;
  124.   TVector3 edge[4]; // vektorji stranic
  125.   double angle_r[4]; // koti ob posameznem vogalu
  126. };
  127. //=================================================================================
  128. // ravnina - krog
  129. class CPlaneR
  130. {
  131.  
  132. public:
  133.   CPlaneR(TVector3 c, TVector3 n0, double R0)
  134.   {center = c; n = n0; _r = R0;};
  135.  
  136.   void Set(TVector3 c, TVector3 n0, double R0)
  137.   {center = c; n = n0; _r = R0;};
  138.  
  139.   int TestIntersection(TVector3 *vec, CRay in);
  140.  
  141.   void Draw(int color = 1, int NN = 32);
  142.  
  143. private:
  144.   TVector3 n, center;
  145.   double _r;
  146. };
  147.  
  148. //=================================================================================
  149. // ravna opticna povrsina: refractor, zrcalo ali povrsina s totalnim odbojem
  150. #define SURF_DUMMY 0
  151. #define SURF_REFRA 1
  152. #define SURF_REFLE 2
  153. #define SURF_TOTAL 3
  154. //#define SURF_TOTAL 1
  155. #define SURF_IMPER 4
  156. //#define SURF_DIFUS 5
  157.  
  158. class CSurface : public CPlane4
  159. {
  160.  
  161. public:
  162.   CSurface(int type0 = 0);
  163.   CSurface(int type0, TVector3 r1, TVector3 r2, TVector3 r3, TVector3 r4,
  164.       double n10, double n20, double reflectivity);
  165.   CSurface(int type0, TVector3 *vr,  double n10, double n20, double reflectivity);
  166.  
  167.   void SetV(TVector3 *vr){CPlane4::Set(vr);};
  168.   void SetType(int type0){type = type0;};
  169.   void SetIndex(double n10, double n20);
  170.   void SetReflection(double reflectivity){reflection = reflectivity;};
  171.   void SetFresnel(int f = 1) {fresnel = f;};
  172.  
  173.   void Set(int type0, TVector3 *vr, double n10, double n20, double reflectivity)
  174.   {type = type0; CPlane4::Set(vr); SetIndex(n10, n20); reflection = reflectivity;};
  175.   void Set(int type0, TVector3 r1, TVector3 r2, TVector3 r3, TVector3 r4, double n10, double n20, double reflectivity)
  176.   {type = type0; CPlane4::Set(r1, r2, r3, r4); SetIndex(n10, n20); reflection = reflectivity;};
  177.  
  178.   int PropagateRay(CRay in, CRay *out, TVector3 *intersection);
  179.  
  180.   double N1_N2(int sign) { return ((sign > 0) ? n1/n2 : n2/n1);};
  181.  
  182.   void Print(){printf("Type = %d\n", type); CPlane4::Print();};
  183.  
  184. private:
  185.   int type; //0 = dummy; 1 = refractor; 2 = reflector; 3 = total reflection
  186.   double n1, n2, n1_n2; //index of refraction, n1 @ +normal, n2 @ -normal
  187.   double reflection; // odbojnost stranic
  188.   TRandom rand; // za racunanje verjetnosti odboja od zrcala
  189.   double cosTtotal; //cosinus mejnega kota totalnega odboja za dana n1 in n2
  190.  
  191.   int fresnel; // ali naj uposteva Fresnelove enacbe
  192. };
  193.  
  194.  
  195. //=================================================================================
  196. //    +-------   n1                
  197. //    |       -------              
  198. //    |          n2  ----+        
  199. //    |                  |        
  200. //   0| M*SiPM           | SiPM    
  201. //    |                  |        
  202. //    |              ----+        
  203. //    |       ------R              
  204. //    +-------                    
  205. //    <-------- d ------->        
  206. //                                
  207.  
  208.  
  209.  
  210. int RCol(int col, int solid);
  211.  
  212.  
  213. class DetectorParameters
  214. {
  215. public:
  216.   DetectorParameters(double a, double b, double d, double active, double n1, double n2, double n3, TVector3 gap, bool coupling):
  217.     _a(a),
  218.     _b(b),
  219.     _d(d),
  220.     _active(active),
  221.     _n1(n1),
  222.     _n2(n2),
  223.     _n3(n3),
  224.     _gap(gap),
  225.     _fresnel(1),
  226.     _guideOn(1),
  227.     offsetY(b/2.0),
  228.     offsetZ(b/2.0),
  229.     _plateOn(1),
  230.     _plateWidth(1),
  231.     _glassOn(0),
  232.     _glassD(0),
  233.     _coupling(coupling)
  234.   {};
  235.   DetectorParameters(double a, double b, double d):
  236.     _a(a),
  237.     _b(b),
  238.     _d(d),
  239.     _active(a),
  240.     _n1(1.0),
  241.     _n2(1.53),
  242.     _n3(1.46),
  243.     _gap(TVector3(0,0,0)),
  244.     _fresnel(1),
  245.     _guideOn(1),
  246.     offsetY(b/2.0),
  247.     offsetZ(b/2.0),
  248.     _plateOn(1),
  249.     _plateWidth(1),
  250.     _glassOn(0),
  251.     _glassD(0),
  252.     _coupling(false)
  253.   {};
  254.   ~DetectorParameters() {};
  255.  
  256.   void setGuide(double a, double b, double d) {
  257.     _a = a;
  258.     _b = b;
  259.     //_M = b/a;
  260.     _d = d;
  261.   };
  262.   void setGap(double x, double y, double z) { _gap = TVector3(x,y,z); };
  263.   void setFresnel(int fresnel) { _fresnel = fresnel; };
  264.   void setGlass(int glassOn, double glassD) { _glassOn = glassOn; _glassD = glassD; };
  265.   void setGuideOn(int guideOn) { _guideOn = guideOn; };
  266.   void setPlate(int plateOn, double plateWidth) { _plateOn = plateOn; _plateWidth = plateWidth; };
  267.   void setIndices(double n1, double n2, double n3) { _n1 = n1; _n2 = n2; _n3 = n3; };
  268.   void setCoupling(bool coupling) { _coupling = coupling; };
  269.  
  270.   double getLightYield() { return ( Power(_b,2)/Power(_active,2)); };
  271.   double getM() { return (_b/_a); };
  272.   double getA() {return _a;};
  273.   double getB() {return _b;};
  274.   double getD() {return _d;};
  275.   double getN1() {return _n1; };
  276.   double getN2() {return _n2;};
  277.   double getN3() {return _n3;};
  278.   double getActive() {return _active;};
  279.   TVector3 getGap() {return _gap;};
  280.   double getPlateWidth() {return _plateWidth; };
  281.   int getFresnel() { return _fresnel; };
  282.   int getGlassOn() { return _glassOn; };
  283.   double getGlassD() { return _glassD; };
  284.   int getGuideOn() { return _guideOn; };
  285.   int getPlateOn() { return _plateOn; };
  286.   double getOffsetY() { return offsetY; };
  287.   double getOffsetZ() { return offsetZ; };
  288.   // bad coupling in the case the small amount of grease was applied
  289.   bool badCoupling() { return _coupling; };
  290.  
  291. private:
  292.   double _a;
  293.   double _b;
  294.   double _d;
  295.   double _active;
  296.   double _n1;
  297.   double _n2;
  298.   double _n3;
  299.   TVector3 _gap;
  300.  
  301.   int _fresnel;
  302.  
  303.   int _guideOn;
  304.   double offsetY;
  305.   double offsetZ;
  306.  
  307.   int _plateOn;
  308.   double _plateWidth;
  309.  
  310.   int _glassOn;
  311.   double _glassD;
  312.   bool _coupling;
  313.  
  314. };
  315.  
  316. class Guide
  317. {
  318. public:
  319.   Guide(TVector3 center0, DetectorParameters& parameters);
  320.   ~Guide() {
  321.     for (int jk=0; jk<6; jk++) delete s_side[jk];
  322.     delete grease;
  323.     delete noCoupling;
  324.     delete hfate;
  325.     delete hnodb_all;
  326.     delete hnodb_exit;
  327.     delete hin;
  328.     delete hout;
  329.   }
  330.  
  331.   Fate PropagateRay(CRay in, CRay *out, int *n_points, TVector3 *points);
  332.  
  333.   double getD() { return _d; }
  334.   double getN1() { return _n1; };
  335.   double getN2() { return _n2; };
  336.   double getN3() { return _n3; };
  337.   TH1F* GetHFate() const {return hfate;};
  338.   TH1F* GetHNOdb_all() const {return hnodb_all;};
  339.   TH1F* GetHNOdb_exit() const {return hnodb_exit;};
  340.   TH2F* GetHIn() const {return hin;};
  341.   TH2F* GetHOut() const {return hout;};
  342.   int GetExitHits() {return (int)hfate->GetBinContent(5);};
  343.   int GetEnteranceHits() {return (int)hfate->GetBinContent(6);};
  344.   void GetVFate(int *out);
  345.   int GetMAXODB() const {return MAX_REFLECTIONS;};
  346.  
  347.   void Draw(int color = 1, int width = 1);
  348.   void DrawSkel(int color = 1, int width = 1);
  349.  
  350. private:
  351.   Fate fate;
  352.   double _d; // parameters
  353.   double _n1, _n2, _n3; // refractive index: n1 above entry surface, n2 inside, n3 after exit
  354.   double _r;
  355.   double _absorption;
  356.   double _A;
  357.   bool _badCoupling;
  358.   TRandom rand; // for material absorption
  359.   CSurface *s_side[6];
  360.   CPlaneR* grease;
  361.   CSurface* noCoupling;
  362.   TVector3 center, vodnik_edge[8];
  363.  
  364.   TH1F *hfate, *hnodb_all, *hnodb_exit;
  365.   TH2F *hin, *hout;
  366. };
  367.  
  368.  
  369.  
  370. class Plate
  371. {
  372. public:
  373.   Plate(DetectorParameters& parameters);
  374.   ~Plate() {
  375.     for (int jk=0; jk<6; jk++) delete sides[jk]; // the same, needs solution
  376.   };
  377.  
  378.   void draw(int color, int width);
  379.   void drawSkel(int color, int width);
  380.   Fate propagateRay(CRay, CRay*, int*, TVector3*);
  381.  
  382. private:
  383.   TVector3 plate_edge[8];
  384.   CSurface *sides[6];
  385. };
  386.  
  387.  
  388. class CDetector
  389. {
  390. public:
  391.   CDetector(TVector3 center0, DetectorParameters& parameters);
  392.   ~CDetector() {
  393.     delete glass;
  394.     delete glass_circle;
  395.     delete hglass;
  396.     delete active;
  397.     delete hactive;
  398.     delete hlaser;
  399.     delete detector;
  400.     delete hdetector;
  401.     delete guide;
  402.     delete plate;
  403.     delete histoPlate;
  404.     //delete window;
  405.     //delete window_circle;
  406.     //delete hwindow;
  407.   }
  408.  
  409.   //void SetLGType(int in = SURF_REFRA, int side = SURF_REFRA, int out = SURF_REFRA)
  410.   //    {type_in = in; type_side = side; type_out = out;};
  411.  
  412.   /*    void SetLG(double SiPM0, double M0, double d0, double n10, double n20, double n30, double R0)
  413.         { SiPM=SiPM0;
  414.                 M=M0;
  415.                 d=d0;
  416.                 n1=n10;
  417.                 n2=n20;
  418.                 n3=n30;
  419.                 R=R0;
  420.                 }; */
  421.   //void SetR(double R0) {R = R0;};
  422.   //void SetGlass(int glass_on0, double glass_d0)
  423.   //{glass_on = glass_on0; glass_d = glass_d0;};
  424.   //void SetGap(double x_gap0, double y_gap0, double z_gap0)
  425.   //{x_gap = x_gap0; y_gap = y_gap0; z_gap = z_gap0;};
  426.   void SetRCol(int in, int lg, int out, int gla)
  427.   {col_in = in; col_lg = lg; col_out = out; col_rgla = gla;};
  428.   void SetDCol(int LG0, int glass0, int active0)
  429.   {col_LG = LG0; col_glass = glass0; col_active = active0;};
  430.   //void SetFresnel(int b)
  431.   //{fresnel = b;};
  432.   //void SetAbsorption(int b, double A0)
  433.   //{absorption = b; A = A0;};
  434.   //void SetWindow(double wR, double d0) {window_R = wR; window_d = d0;};
  435.   //void SetGuideOn(int b)
  436.   //{guide_on = b;};
  437.  
  438.   //void Init();
  439.  
  440.   int Propagate(CRay, CRay*, int);
  441.  
  442.   Guide* GetLG() const {return guide;};
  443.   //TH2F* GetHWindow() const {return hwindow;};
  444.   TH2F* GetHGlass() const {return hglass;};
  445.   TH2F* GetHActive() const {return hactive;};
  446.   TH2F* GetHLaser() const {return hlaser;};
  447.   TH2F* GetHDetector() const {return hdetector;};
  448.   TH2F* GetHPlate() const {return histoPlate;};
  449.  
  450.   //double GetSiPM() const {return SiPM;};
  451.   //double GetM() const {return M;};
  452.   //double GetD() const {return d;};
  453.   //double GetR() const {return R;};
  454.   //TVector3 GetVGap() const {return (TVector3(x_gap, y_gap, z_gap));};
  455.   //int IsGlass() {return glass_on;};
  456.   //double getGlassWidth() {return glass_d;};
  457.   //double GetWindowR() {return window_R;};
  458.   //double GetWindowD() {return window_d;};
  459.  
  460.   void Draw(int width = 2);
  461.  
  462. private:
  463.   Fate fate;
  464.   TVector3 center;
  465.  
  466.   //int type_in, type_side, type_out;
  467.   //double SiPM, M, d, n1, n2, n3, R;
  468.   //double detectorActive;
  469.  
  470.  
  471.   int glass_on;
  472.   double glass_d;
  473.   CSurface *glass;
  474.   CPlaneR *glass_circle;
  475.   TH2F *hglass;
  476.  
  477.   //double x_gap, y_gap, z_gap;
  478.   CPlane4* active;
  479.   CPlaneR* grease;
  480.   TH2F *hactive, *hlaser;
  481.  
  482.   CPlane4 *detector;
  483.   TH2F *hdetector;
  484.  
  485.   int col_in, col_lg, col_out, col_rgla;
  486.   int col_LG, col_glass, col_active;
  487.  
  488.   //int fresnel, absorption;
  489.   //double A;
  490.  
  491.   int guide_on;
  492.  
  493.   //double window_R, window_d;
  494.   //CSurface *window;
  495.   //CPlaneR *window_circle;
  496.   //TH2F *hwindow;
  497.  
  498.   Guide *guide;
  499.   Plate *plate;
  500.  
  501.   double _plateWidth;
  502.   int _plateOn;
  503.  
  504.   TH2F *histoPlate;
  505.  
  506.   double offsetY;
  507.   double offsetZ;
  508. };
  509.  
  510.  
  511.  
  512. #endif
  513.