Subversion Repositories f9daq

Rev

Rev 54 | Go to most recent revision | Details | Last modification | View Log | RSS feed

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