Subversion Repositories f9daq

Compare Revisions

Ignore whitespace Rev 70 → Rev 71

/lightguide/trunk/include/guide.h
441,7 → 441,7
//double GetR() const {return R;};
//TVector3 GetVGap() const {return (TVector3(x_gap, y_gap, z_gap));};
//int IsGlass() {return glass_on;};
//double GetGlass() {return glass_d;};
//double getGlassWidth() {return glass_d;};
//double GetWindowR() {return window_R;};
//double GetWindowD() {return window_d;};
 
/lightguide/trunk/src/raySimulator.cpp
97,9 → 97,9
int cpc = 1;
cdata->cd(cpc++); (detector->histoPlate)->Draw("COLZ");
cdata->cd(cpc++); ((detector->GetLG())->GetHIn())->Draw("COLZ");
cdata->cd(cpc++); ((detector->GetLG())->GetHOut())->Draw("COLZ");
cdata->cd(cpc++); //((detector->GetLG())->GetHOut())->Draw("COLZ");
(detector->GetHGlass())->Draw("COLZ");
cdata->cd(cpc++); (detector->GetHActive())->Draw("COLZ");
cdata->cd(cpc++); (detector->GetHLaser())->Draw("COLZ");
cdata->cd(cpc++); (detector->GetHDetector())->Draw("COLZ");
/lightguide/trunk/src/guide.cpp
385,7 → 385,8
int CSurface::PropagateRay(CRay in, CRay *out, TVector3 *intersection)
{
if (dbg) printf("--- CSurface::PropagateRay ---\n");
double cosTi, cosTt, p_ref;
double cosTi; // incident ray angle
double cosTt; // transmited ray angle
TVector3 intersect, transmit;
 
if( !(GetIntersection(&intersect, in) == 1) )
437,6 → 438,9
}
// ----------------------------------------------------------------------------
// reflection probability
double p_ref = rand.Uniform(0.0, 1.0);
if(type == SURF_TOTAL) type = SURF_REFRA;
switch(type){
// ----------------------------------------------------------------------------
458,11 → 462,10
if(dbg) printf(" cosTtotal = %lf (Ttotal = %lf)\n", cosTtotal, TMath::ACos(cosTtotal)*DEGREE);
// reflection dependance on polarization missing
// reflection hardcoded to 0.96
p_ref = rand.Uniform(0.0, 1.0);
if (dbg) printf(" reflection probability = %f\n", p_ref);
// If n1>n2 and theta>thetaCritical, total reflection
if( (cosTi < cosTtotal) && (p_ref < reflection) ) { // totalni odboj z verjetnostjo "reflection"
if(cosTi < cosTtotal) {
if(dbg) printf(" TOTAL\n");
transmit = in.GetN() + sign_n*2*cosTi*n;
484,38 → 487,38
transmit = N1_N2(sign_n)*in.GetN() + sign_n*(N1_N2(sign_n)*cosTi - cosTt)*n;
if(dbg) {printf(" transmit.Unit() = "); printv(transmit.Unit());}
 
if(dbg) {
cosTN = transmit.Unit() * n;
printf(" cosTN = %lf (TN = %lf) (Abs(TN) = %lf)\n", cosTN, TMath::ACos(cosTN)*DEGREE, TMath::ACos(TMath::Abs(cosTN))*DEGREE);
}
}
//if(cosTi<=cosTtotal) cosTt = TMath::Sqrt(1 - TMath::Power(N1_N2(sign_n), 2)*(1 - TMath::Power(cosTi, 2)));
//if(fresnel) {
r_te = (n1*cosTi - n2*cosTt)/(n1*cosTi + n2*cosTt); // transverse
r_tm = (n2*cosTi - n1*cosTt)/(n1*cosTt + n2*cosTi); // paralel
r_te = (n1*cosTi - n2*cosTt)/(n1*cosTi + n2*cosTt); // transverse
r_tm = (n2*cosTi - n1*cosTt)/(n1*cosTt + n2*cosTi); // paralel
if(dbg) printf(" r_te = %lf, r_tm = %lf\n", r_te, r_tm);
if(dbg) printf(" r_te = %lf, r_tm = %lf\n", r_te, r_tm);
// transmited polarization
v_tm_t = -v_te.Cross(transmit);
v_tm_t = v_tm_t.Unit();
pol_t = a_te * (1.0 - TMath::Abs(r_te)) * v_te + a_tm * (1.0 - TMath::Abs(r_tm)) * v_tm_t;
// transmited polarization
v_tm_t = -v_te.Cross(transmit);
v_tm_t = v_tm_t.Unit();
pol_t = a_te * (1.0 - TMath::Abs(r_te)) * v_te + a_tm * (1.0 - TMath::Abs(r_tm)) * v_tm_t;
if(dbg) {
if(dbg) {
printf(" v_tm_t = "); printv(v_tm_t);
printf(" pol_t = "); printv(pol_t);
}
}
// Fresnel coefficients
R_te = TMath::Power(r_te, 2);
R_tm = TMath::Power(r_tm, 2);
R_f = a_te*a_te*R_te + a_tm*a_tm*R_tm;
// Fresnel coefficients
R_te = TMath::Power(r_te, 2);
R_tm = TMath::Power(r_tm, 2);
R_f = a_te*a_te*R_te + a_tm*a_tm*R_tm;
if (dbg) printf(" R_te = %lf, R_tm = %lf, R_f = %lf\n", R_te, R_tm, R_f);
}
if (dbg) printf(" R_te = %lf, R_tm = %lf, R_f = %lf\n", R_te, R_tm, R_f);
}
if(p_ref >= R_f) { // se lomi
if (dbg) printf(" SURFACE REFRACTED. Return.\n");
if(p_ref >= R_f) { // se lomi
if (dbg) printf(" SURFACE REFRACTED. Return.\n");
out->Set(intersect, transmit);
out->SetPolarization(pol_t);
return REFRACTION;
919,7 → 922,7
// additional glass between at top of SiPM
// example: epoxy n=1.60
double n4 = 1.60;
double n4 = 1.57;
TVector3 plane_v[4];
int nBins = nch + 1;
double p_size = b/2.0;
1127,7 → 1130,8
if (dbg) printf("Detector: fate = hit or refracted");
ray1 = ray0;
if(draw) {
if(glass_on) ray1->Draw(center.x(), center.x() /*+ window_d*/);
//double epoxy = parameters->getGlassD();
if(glass_on) ray1->Draw(center.x(), center.x() + glass_d);
else ray1->DrawS(center.x(), 10.0);
}
}
1134,22 → 1138,25
fate = missed; // zgresil aktivno povrsino
if(glass_on) {
*ray0 = *ray1; ray1->SetColor(col_rgla);
*ray0 = *ray1;
ray1->SetColor(col_rgla);
fate_glass = glass->PropagateRay(*ray0, ray1, &presecisce);
if(fate_glass == 1) {
if(fate_glass == REFRACTION) {
hglass->Fill(presecisce.y(), presecisce.z());
if(draw) ray1->DrawS(presecisce.x(), 10.0);
if(active->TestIntersection(&presecisce, *ray1)) {
fate = hitExit;
hactive->Fill(offsetY + presecisce.y(), offsetZ + presecisce.z());
hlaser->Fill((in.GetR()).y(), (in.GetR()).z());
}
if(detector->TestIntersection(&presecisce, *ray1))
hdetector->Fill(presecisce.y(), presecisce.z());
} else if(fate_glass == 2) {
//if(active->TestIntersection(&presecisce, *ray1)) {
//fate = hitExit;
//hactive->Fill(offsetY + presecisce.y(), offsetZ + presecisce.z());
//hlaser->Fill((in.GetR()).y() + offsetY, (in.GetR()).z() + offsetZ);
//}
//if(detector->TestIntersection(&presecisce, *ray1))
//hdetector->Fill(offsetY + presecisce.y(), offsetZ + presecisce.z());
//} else if(fate_glass == REFLECTION) {
else
if(draw) ray1->DrawS(presecisce.x(), 10.0);
}
} else {
}
 
// Main test: ray and SiPM surface
if(active->TestIntersection(&presecisce, *ray1)) {
fate = hitExit;
1159,7 → 1166,7
// If it is on the same plane as SiPM
if(detector->TestIntersection(&presecisce, *ray1))
hdetector->Fill(offsetY + presecisce.y(), offsetZ + presecisce.z());
}
//}
//} else {
//if(draw) ray1->Draw(presecisce.x(), center.x()+d+window_d);
//}