#include "TGeoMaterial.h"
#include "TGeoMedium.h"
#include "TGeoManager.h"
#include "TGeoCompositeShape.h"
#include "TGeoCone.h"
#include "TGeoMatrix.h"
#include "TGeoBoolNode.h"
#include "TGeoXtru.h"
#include "TRandom.h"
#include "TPolyLine3D.h"
#include "Riostream.h"

#define HEX_LADO 354.49/2.0
#define HEX_ALTU 307.00/2.0
#define HEX_VERT       66.4

#define VIDRO_NHSEG    50
#define VIDRO_NVSEG    20
#define PMT_RAIO    103.5
#define PMT_HEIGHT   82.9

#define COS30 cos(3.141592654/6.0)
#define STEP_Y HEX_ALTU*2.0
#define STEP_X STEP_Y*COS30

#define DET_W1  1152.10
#define DET_W2  1381.5
#define DET_H   (900.0 + 28.46)
#define OFF_X   (DET_W1/2.0 - HEX_LADO)
#define OFF_Y1  (DET_W2/2.0 - HEX_ALTU)
#define OFF_Y2  (DET_W2/2.0 -   STEP_Y)

void Reflexao()
{
	// cria mundo ---------------------------------------------------------

	TGeoMaterial *matagua = new TGeoMaterial("Agua"   ,0,0,0);
	TGeoMaterial *matmirr = new TGeoMaterial("Espelho",0,0,0);
	TGeoMaterial *matglas = new TGeoMaterial("Vidro"  ,0,0,0);

	TGeoMedium   *medagua = new TGeoMedium("Agua"   ,1, matagua);
	TGeoMedium   *medmirr = new TGeoMedium("Espelho",1, matmirr);
	TGeoMedium   *medglas = new TGeoMedium("Vidro"  ,1, matglas);

	// caixa do detector --------------------------------------------------

	TGeoVolume *det = gGeoManager->MakeBox("Detector", medagua, DET_W1/2.0, DET_W2/2.0 + HEX_VERT, DET_H/2.0 + HEX_VERT);

	// espelho inferior ---------------------------------------------------

	TGeoShape *sbox  = new TGeoBBox("mbox", DET_W1/2.0, DET_W2/2.0, HEX_VERT/2.0);
	TGeoShape *scone = new TGeoCone("cone", HEX_VERT/2.0*1.2, 0.0, PMT_RAIO, 0.0, HEX_LADO);
	
	TGeoTranslation    *t1 = new TGeoTranslation("t1", -OFF_X, -OFF_Y1, 0);
	TGeoSubtraction    *s1 = new TGeoSubtraction(sbox, scone, NULL, t1);
	TGeoCompositeShape *c1 = new TGeoCompositeShape("c1", s1);

	TGeoTranslation    *t2 = new TGeoTranslation("t2", -OFF_X, -OFF_Y1+STEP_Y, 0);
	TGeoSubtraction    *s2 = new TGeoSubtraction(c1, scone, NULL, t2);
	TGeoCompositeShape *c2 = new TGeoCompositeShape("c2", s2);

	TGeoTranslation    *t3 = new TGeoTranslation("t3", -OFF_X, -OFF_Y1+2.0*STEP_Y, 0);
	TGeoSubtraction    *s3 = new TGeoSubtraction(c2, scone, NULL, t3);
	TGeoCompositeShape *c3 = new TGeoCompositeShape("c3", s3);

	TGeoTranslation    *t4 = new TGeoTranslation("t4", -OFF_X, -OFF_Y1+3.0*STEP_Y, 0);
	TGeoSubtraction    *s4 = new TGeoSubtraction(c3, scone, NULL, t4);
	TGeoCompositeShape *c4 = new TGeoCompositeShape("c4", s4);

	TGeoTranslation    *t5 = new TGeoTranslation("t5", -OFF_X+STEP_X, -OFF_Y2, 0);
	TGeoSubtraction    *s5 = new TGeoSubtraction(c4, scone, NULL, t5);
	TGeoCompositeShape *c5 = new TGeoCompositeShape("c5", s5);

	TGeoTranslation    *t6 = new TGeoTranslation("t6", -OFF_X+STEP_X, -OFF_Y2+STEP_Y, 0);
	TGeoSubtraction    *s6 = new TGeoSubtraction(c5, scone, NULL, t6);
	TGeoCompositeShape *c6 = new TGeoCompositeShape("c6", s6);

	TGeoTranslation    *t7 = new TGeoTranslation("t7", -OFF_X+STEP_X, -OFF_Y2+2.0*STEP_Y, 0);
	TGeoSubtraction    *s7 = new TGeoSubtraction(c6, scone, NULL, t7);
	TGeoCompositeShape *c7 = new TGeoCompositeShape("c7", s7);

	TGeoTranslation    *t8 = new TGeoTranslation("t8", -OFF_X+STEP_X, -OFF_Y2+3.0*STEP_Y, 0);
	TGeoSubtraction    *s8 = new TGeoSubtraction(c7, scone, NULL, t8);
	TGeoCompositeShape *c8 = new TGeoCompositeShape("c8", s8);

	TGeoTranslation    *t9 = new TGeoTranslation("t9", -OFF_X+2.0*STEP_X, -OFF_Y1, 0);
	TGeoSubtraction    *s9 = new TGeoSubtraction(c8, scone, NULL, t9);
	TGeoCompositeShape *c9 = new TGeoCompositeShape("c9", s9);

	TGeoTranslation    *t10 = new TGeoTranslation("t10", -OFF_X+2.0*STEP_X, -OFF_Y1+STEP_Y, 0);
	TGeoSubtraction    *s10 = new TGeoSubtraction(c9, scone, NULL, t10);
	TGeoCompositeShape *c10 = new TGeoCompositeShape("c10", s10);

	TGeoTranslation    *t11 = new TGeoTranslation("t11", -OFF_X+2.0*STEP_X, -OFF_Y1+2.0*STEP_Y, 0);
	TGeoSubtraction    *s11 = new TGeoSubtraction(c10, scone, NULL, t11);
	TGeoCompositeShape *c11 = new TGeoCompositeShape("c11", s11);

	TGeoTranslation    *t12 = new TGeoTranslation("t12", -OFF_X+2.0*STEP_X, -OFF_Y1+3.0*STEP_Y, 0);
	TGeoSubtraction    *s12 = new TGeoSubtraction(c11, scone, NULL, t12);
	TGeoCompositeShape *c12 = new TGeoCompositeShape("c12", s12);

	TGeoTranslation    *t13 = new TGeoTranslation("t13", -OFF_X+3.0*STEP_X, -OFF_Y2, 0);
	TGeoSubtraction    *s13 = new TGeoSubtraction(c12, scone, NULL, t13);
	TGeoCompositeShape *c13 = new TGeoCompositeShape("c13", s13);

	TGeoTranslation    *t14 = new TGeoTranslation("t14", -OFF_X+3.0*STEP_X, -OFF_Y2+STEP_Y, 0);
	TGeoSubtraction    *s14 = new TGeoSubtraction(c13, scone, NULL, t14);
	TGeoCompositeShape *c14 = new TGeoCompositeShape("c14", s14);

	TGeoTranslation    *t15 = new TGeoTranslation("t15", -OFF_X+3.0*STEP_X, -OFF_Y2+2.0*STEP_Y, 0);
	TGeoSubtraction    *s15 = new TGeoSubtraction(c14, scone, NULL, t15);
	TGeoCompositeShape *c15 = new TGeoCompositeShape("c15", s15);

	TGeoTranslation    *t16 = new TGeoTranslation("t16", -OFF_X+3.0*STEP_X, -OFF_Y2+3.0*STEP_Y, 0);
	TGeoSubtraction    *s16 = new TGeoSubtraction(c15, scone, NULL, t16);
	TGeoCompositeShape *c16 = new TGeoCompositeShape("c16", s16);
	TGeoVolume         *v16 = new TGeoVolume        ("v16", c16, medmirr);

	v16->SetLineColor(kGray-1);

	TGeoTranslation *t17 = new TGeoTranslation("t17", 0.0,   0.0, -(DET_H/2.0 + HEX_VERT/2.0));
	TGeoRotation     *r1 = new TGeoRotation   ("r1" , 0.0, 180.0, 0.0);
	TGeoCombiTrans  *ct1 = new TGeoCombiTrans ("ct1", 0.0,   0.0,   DET_H/2.0 + HEX_VERT/2.0, r1);
	det->AddNode(v16, 1, t17);
	det->AddNode(v16, 2, ct1);

	// espelho lateral ----------------------------------------------------

	TGeoBBox *slbox = new TGeoBBox("lbox", DET_W1/2.0, DET_H/2.0 + HEX_VERT, HEX_VERT/2.0);

	TGeoTranslation    *t18 = new TGeoTranslation("t18", HEX_ALTU, 0, 0);
	TGeoSubtraction    *s18 = new TGeoSubtraction(slbox, scone, NULL, t18);
	TGeoCompositeShape *c18 = new TGeoCompositeShape("c18", s18);

	TGeoTranslation    *t19 = new TGeoTranslation("t19", -HEX_ALTU, 0, 0);
	TGeoSubtraction    *s19 = new TGeoSubtraction(c18, scone, NULL, t19);
	TGeoCompositeShape *c19 = new TGeoCompositeShape("c19", s19);

	TGeoTranslation    *t20 = new TGeoTranslation("t20", 0, STEP_X, 0);
	TGeoSubtraction    *s20 = new TGeoSubtraction(c19, scone, NULL, t20);
	TGeoCompositeShape *c20 = new TGeoCompositeShape("c20", s20);

	TGeoTranslation    *t21 = new TGeoTranslation("t21", 0, -STEP_X, 0);
	TGeoSubtraction    *s21 = new TGeoSubtraction(c20, scone, NULL, t21);
	TGeoCompositeShape *c21 = new TGeoCompositeShape("c21", s21);
	TGeoVolume         *v21 = new TGeoVolume        ("v21", c21, medmirr);

	v21->SetLineColor(kGray-1);

	TGeoRotation    *r22 = new TGeoRotation("r22", 0.0, 90.0, 0.0);
	TGeoCombiTrans *ct22 = new TGeoCombiTrans("ct22", 0.0, DET_W2/2.0 + HEX_VERT/2.0, 0.0, r22);
	det->AddNode(v21, 1, ct22);

	TGeoTranslation    *t23 = new TGeoTranslation("t23", STEP_X, 0, 0);
	TGeoSubtraction    *s23 = new TGeoSubtraction(slbox, scone, NULL, t23);
	TGeoCompositeShape *c23 = new TGeoCompositeShape("c23", s23);

	TGeoTranslation    *t24 = new TGeoTranslation("t24", -STEP_X, 0, 0);
	TGeoSubtraction    *s24 = new TGeoSubtraction(c23, scone, NULL, t24);
	TGeoCompositeShape *c24 = new TGeoCompositeShape("c24", s24);

	TGeoTranslation    *t25 = new TGeoTranslation("t25", 0, HEX_ALTU, 0);
	TGeoSubtraction    *s25 = new TGeoSubtraction(c24, scone, NULL, t25);
	TGeoCompositeShape *c25 = new TGeoCompositeShape("c25", s25);

	TGeoTranslation    *t26 = new TGeoTranslation("t26", 0, -HEX_ALTU, 0);
	TGeoSubtraction    *s26 = new TGeoSubtraction(c25, scone, NULL, t26);
	TGeoCompositeShape *c26 = new TGeoCompositeShape("c26", s26);
	TGeoVolume         *v26 = new TGeoVolume        ("v26", c26, medmirr);

	v26->SetLineColor(kGray-1);

	TGeoRotation   *r27   = new TGeoRotation("r27", 0.0, -90.0, 0.0);
	TGeoCombiTrans *ctpmt = new TGeoCombiTrans("ctpmt", 0.0, -DET_W2/2.0 - HEX_VERT/2.0, 0.0, r27);
	det->AddNode(v26, 1, ctpmt);

	// cria PMT -----------------------------------------------------------

	TGeoVolume *Vidro = gGeoManager->MakeXtru("Vidro", medglas, VIDRO_NVSEG);
	TGeoXtru    *xtru = (TGeoXtru*)Vidro->GetShape();

	Int_t i;
	Double_t x[VIDRO_NHSEG], y[VIDRO_NHSEG];
	Double_t dtheta = 2.0*3.1415926/VIDRO_NHSEG;
	Double_t  theta = 0.0;
	
	for (i = 0; i < VIDRO_NHSEG; i++)
	{
		x[i] = PMT_RAIO*sin(theta);
		y[i] = PMT_RAIO*cos(theta);
		theta += dtheta;
	}

	xtru->DefinePolygon(VIDRO_NHSEG,x,y);

	Double_t h  = -HEX_VERT/2.0;
	Double_t ds = 1.0/VIDRO_NVSEG;
	Double_t s  = 1.0;
	for (int i = 0; i < VIDRO_NVSEG; i++)
	{
		xtru->DefineSection(i, h, 0.0, 0.0, s);
		s -= ds;  
		h = sqrt(PMT_HEIGHT*PMT_HEIGHT*(1.0 - s*s)) -HEX_VERT/2.0;
	}

	Vidro->SetTransparency(30);
	Vidro->SetLineColor(kOrange);

	// posiciona PMTs -----------------------------------------------------

	TGeoTranslation    *tpmt1  = new TGeoTranslation("tpmt", -OFF_X,            -OFF_Y1,            -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro,  1, tpmt1);
	TGeoTranslation    *tpmt2  = new TGeoTranslation("tpmt", -OFF_X,            -OFF_Y1+    STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro,  2, tpmt2);
	TGeoTranslation    *tpmt3  = new TGeoTranslation("tpmt", -OFF_X,            -OFF_Y1+2.0*STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro,  3, tpmt3);
	TGeoTranslation    *tpmt4  = new TGeoTranslation("tpmt", -OFF_X,            -OFF_Y1+3.0*STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro,  4, tpmt4);
	TGeoTranslation    *tpmt5  = new TGeoTranslation("tpmt", -OFF_X+    STEP_X, -OFF_Y2,            -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro,  5, tpmt5);
	TGeoTranslation    *tpmt6  = new TGeoTranslation("tpmt", -OFF_X+    STEP_X, -OFF_Y2+    STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro,  6, tpmt6);
	TGeoTranslation    *tpmt7  = new TGeoTranslation("tpmt", -OFF_X+    STEP_X, -OFF_Y2+2.0*STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro,  7, tpmt7);
	TGeoTranslation    *tpmt8  = new TGeoTranslation("tpmt", -OFF_X+    STEP_X, -OFF_Y2+3.0*STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro,  8, tpmt8);
	TGeoTranslation    *tpmt9  = new TGeoTranslation("tpmt", -OFF_X+2.0*STEP_X, -OFF_Y1,            -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro,  9, tpmt9);
	TGeoTranslation    *tpmt10 = new TGeoTranslation("tpmt", -OFF_X+2.0*STEP_X, -OFF_Y1+    STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro, 10, tpmt10);
	TGeoTranslation    *tpmt11 = new TGeoTranslation("tpmt", -OFF_X+2.0*STEP_X, -OFF_Y1+2.0*STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro, 11, tpmt11);
	TGeoTranslation    *tpmt12 = new TGeoTranslation("tpmt", -OFF_X+2.0*STEP_X, -OFF_Y1+3.0*STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro, 12, tpmt12);
	TGeoTranslation    *tpmt13 = new TGeoTranslation("tpmt", -OFF_X+3.0*STEP_X, -OFF_Y2,            -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro, 13, tpmt13);
	TGeoTranslation    *tpmt14 = new TGeoTranslation("tpmt", -OFF_X+3.0*STEP_X, -OFF_Y2+    STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro, 14, tpmt14);
	TGeoTranslation    *tpmt15 = new TGeoTranslation("tpmt", -OFF_X+3.0*STEP_X, -OFF_Y2+2.0*STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro, 15, tpmt15);
	TGeoTranslation    *tpmt16 = new TGeoTranslation("tpmt", -OFF_X+3.0*STEP_X, -OFF_Y2+3.0*STEP_Y, -(DET_H/2.0 + HEX_VERT/2.0));	det->AddNode(Vidro, 16, tpmt16);

	TGeoRotation   *rpmt1  = new TGeoRotation  ("rpmt", 0.0, 180.0, 0.0);
	TGeoCombiTrans *cpmt1  = new TGeoCombiTrans("cpmt", -OFF_X,            -OFF_Y2,            DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 17, cpmt1);
	TGeoCombiTrans *cpmt2  = new TGeoCombiTrans("cpmt", -OFF_X,            -OFF_Y2+    STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 18, cpmt2);
	TGeoCombiTrans *cpmt3  = new TGeoCombiTrans("cpmt", -OFF_X,            -OFF_Y2+2.0*STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 19, cpmt3);
	TGeoCombiTrans *cpmt4  = new TGeoCombiTrans("cpmt", -OFF_X,            -OFF_Y2+3.0*STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 20, cpmt4);
	TGeoCombiTrans *cpmt5  = new TGeoCombiTrans("cpmt", -OFF_X+    STEP_X, -OFF_Y1,            DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 21, cpmt5);
	TGeoCombiTrans *cpmt6  = new TGeoCombiTrans("cpmt", -OFF_X+    STEP_X, -OFF_Y1+    STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 22, cpmt6);
	TGeoCombiTrans *cpmt7  = new TGeoCombiTrans("cpmt", -OFF_X+    STEP_X, -OFF_Y1+2.0*STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 23, cpmt7);
	TGeoCombiTrans *cpmt8  = new TGeoCombiTrans("cpmt", -OFF_X+    STEP_X, -OFF_Y1+3.0*STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 24, cpmt8);
	TGeoCombiTrans *cpmt9  = new TGeoCombiTrans("cpmt", -OFF_X+2.0*STEP_X, -OFF_Y2,            DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 25, cpmt9);
	TGeoCombiTrans *cpmt10 = new TGeoCombiTrans("cpmt", -OFF_X+2.0*STEP_X, -OFF_Y2+    STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 26, cpmt10);
	TGeoCombiTrans *cpmt11 = new TGeoCombiTrans("cpmt", -OFF_X+2.0*STEP_X, -OFF_Y2+2.0*STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 27, cpmt11);
	TGeoCombiTrans *cpmt12 = new TGeoCombiTrans("cpmt", -OFF_X+2.0*STEP_X, -OFF_Y2+3.0*STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 28, cpmt12);
	TGeoCombiTrans *cpmt13 = new TGeoCombiTrans("cpmt", -OFF_X+3.0*STEP_X, -OFF_Y1,            DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 29, cpmt13);
	TGeoCombiTrans *cpmt14 = new TGeoCombiTrans("cpmt", -OFF_X+3.0*STEP_X, -OFF_Y1+    STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 30, cpmt14);
	TGeoCombiTrans *cpmt15 = new TGeoCombiTrans("cpmt", -OFF_X+3.0*STEP_X, -OFF_Y1+2.0*STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 31, cpmt15);
	TGeoCombiTrans *cpmt16 = new TGeoCombiTrans("cpmt", -OFF_X+3.0*STEP_X, -OFF_Y1+3.0*STEP_Y, DET_H/2.0 + HEX_VERT/2.0, rpmt1); det->AddNode(Vidro, 32, cpmt16);
	
	TGeoRotation   *rpmt2   = new TGeoRotation  ("rpmt", 0.0, -90.0, 0.0);
	TGeoCombiTrans *cpmt17  = new TGeoCombiTrans("cpmt",     0.0, -DET_W2/2.0 - HEX_VERT/2.0,  HEX_ALTU, rpmt2); det->AddNode(Vidro, 33, cpmt17);
	TGeoCombiTrans *cpmt18  = new TGeoCombiTrans("cpmt",     0.0, -DET_W2/2.0 - HEX_VERT/2.0, -HEX_ALTU, rpmt2); det->AddNode(Vidro, 34, cpmt18);
	TGeoCombiTrans *cpmt19  = new TGeoCombiTrans("cpmt",  STEP_X, -DET_W2/2.0 - HEX_VERT/2.0,       0.0, rpmt2); det->AddNode(Vidro, 35, cpmt19);
	TGeoCombiTrans *cpmt20  = new TGeoCombiTrans("cpmt", -STEP_X, -DET_W2/2.0 - HEX_VERT/2.0,       0.0, rpmt2); det->AddNode(Vidro, 36, cpmt20);

	TGeoRotation   *rpmt3   = new TGeoRotation  ("rpmt", 0.0,  90.0, 0.0);
	TGeoCombiTrans *cpmt21  = new TGeoCombiTrans("cpmt",       0.0,  DET_W2/2.0 + HEX_VERT/2.0,  STEP_X, rpmt3); det->AddNode(Vidro, 37, cpmt21);
	TGeoCombiTrans *cpmt22  = new TGeoCombiTrans("cpmt",       0.0,  DET_W2/2.0 + HEX_VERT/2.0, -STEP_X, rpmt3); det->AddNode(Vidro, 38, cpmt22);
	TGeoCombiTrans *cpmt23  = new TGeoCombiTrans("cpmt",  HEX_ALTU,  DET_W2/2.0 + HEX_VERT/2.0,     0.0, rpmt3); det->AddNode(Vidro, 39, cpmt23);
	TGeoCombiTrans *cpmt24  = new TGeoCombiTrans("cpmt", -HEX_ALTU,  DET_W2/2.0 + HEX_VERT/2.0,     0.0, rpmt3); det->AddNode(Vidro, 40, cpmt24);

	// fecha a geometria --------------------------------------------------

	gGeoManager->SetTopVolume(det);
	gGeoManager->CloseGeometry();
	//gGeoManager->SetTopVisible();
	//det->SetTransparency(50);
	det->Draw("ogl");

	// Inicializa uma trajetoria aleatoriamente ---------------------------

	Double_t *pos; pos = new Double_t[3];
	pos[0] = gRandom->Uniform(-400.0, 400.0); gRandom->Rndm();
	pos[1] = gRandom->Uniform(-400.0, 400.0); gRandom->Rndm();
	pos[2] = gRandom->Uniform(-400.0, 400.0); gRandom->Rndm();

	Double_t *dir; dir = new Double_t[3];
	dir[0] = gRandom->Uniform(-1.0, 1.0); gRandom->Rndm();
	dir[1] = gRandom->Uniform(-1.0, 1.0); gRandom->Rndm();
	dir[2] = gRandom->Uniform(-1.0, 1.0); gRandom->Rndm();

	// Reflete photon ate encontrar um vidro de PMT -----------------------

	Double_t  m, cosa;
	Double_t *n;
	Int_t     N = -1;
	Double_t  comp_tot = 0.0;

	bool saiu = kFALSE;

	while (!saiu)
	{
		// normaliza direcao ------------------------------------------

		m = sqrt(dir[0]*dir[0] + dir[1]*dir[1] + dir[2]*dir[2]);
		dir[0] = dir[0]/m; dir[1] = dir[1]/m; dir[2] = dir[2]/m;

		// inicializa track -------------------------------------------

		gGeoManager->InitTrack(pos, dir);

		TPolyLine3D *pl = new TPolyLine3D(2);
		pl->SetLineColor(kRed);
		pl->SetPoint(0, pos[0], pos[1], pos[2]);

		// propaga track ate proximo objeto ---------------------------

		gGeoManager->FindNextBoundaryAndStep();
		n = gGeoManager->FindNormalFast();
		comp_tot += gGeoManager->GetStep();
		N++;

		pos = (Double_t *) gGeoManager->GetCurrentPoint();
		pos[0] -= 1E-6*dir[0];
        	pos[1] -= 1E-6*dir[1];
        	pos[2] -= 1E-6*dir[2];
		if (gGeoManager->GetStep() <= 1e-06) cout << "perdeu" << endl;

		// pega segundo ponto do track --------------------------------

		pl->SetPoint(1, pos[0], pos[1], pos[2]);
		pl->Draw();

		// pega e desenha normal no novo ponto ------------------------

		TPolyLine3D *norm = new TPolyLine3D(2);
		norm->SetLineColor(kGreen);
		norm->SetPoint(0, pos[0], pos[1], pos[2]);
		norm->SetPoint(1, pos[0] - n[0]*50.0, pos[1] - n[1]*50.0, pos[2] - n[2]*50.0);
		norm->Draw();

		// computa proxima direcao de propagacao ----------------------

		cosa = dir[0]*n[0]+dir[1]*n[1]+dir[2]*n[2];
		dir[0] = dir[0] - 2.0*cosa*n[0];
		dir[1] = dir[1] - 2.0*cosa*n[1];
		dir[2] = dir[2] - 2.0*cosa*n[2];

		gGeoManager->cd("/Detector_1");

		if (strcmp(gGeoManager->GetNextNode()->GetVolume()->GetName(), "Vidro") == 0) saiu = kTRUE;
	}

	cout << endl;
	cout << "Numero de reflexoes = " << N << endl;
	cout << "Comprimento total   = " << comp_tot/1000.0 << " metros." << endl;
	cout << endl;
}