8 |
|
|
9 |
|
#include "standard.h" |
10 |
|
|
11 |
+ |
#include "ray.h" |
12 |
+ |
#include "calcomp.h" |
13 |
+ |
#include "func.h" |
14 |
+ |
|
15 |
|
#ifdef NOSTRUCTASSIGN |
16 |
|
static double err = "No structure assignment!"; /* generate compiler error */ |
17 |
|
#endif |
18 |
|
|
19 |
|
|
20 |
|
static double |
21 |
< |
Sqrt(x) |
22 |
< |
double x; |
21 |
> |
Sqrt( |
22 |
> |
double x |
23 |
> |
) |
24 |
|
{ |
25 |
|
if (x < 0.) |
26 |
|
return(0.); |
84 |
|
static double fact_ref[4]={1.0,1.0,1.0,1.0}; /* facteurs de reflexion */ |
85 |
|
static double tolerance; /* degre de tol. pour les amalgames */ |
86 |
|
static double tolsource; /* degre de tol. pour les sources */ |
82 |
– |
static double Nx; |
87 |
|
static int bidon; |
88 |
|
#define BADVAL (-10) |
89 |
|
static long prismclock = -1; |
93 |
|
static int nbrayons; /* indice des rayons sortants */ |
94 |
|
static TRAYON *ray; /* tableau des rayons sortants */ |
95 |
|
static TRAYON *raytemp; /* variable temporaire */ |
92 |
– |
static TRAYON rtemp; /* variable temporaire */ |
96 |
|
|
94 |
– |
extern double argument(); |
95 |
– |
extern double varvalue(); |
96 |
– |
extern double funvalue(); |
97 |
– |
extern long eclock; |
97 |
|
|
98 |
+ |
static void prepare_matrices(void); |
99 |
+ |
static void tfm(MAT4 mat, FVECT v_old, FVECT v_new); |
100 |
+ |
static double prob_alpha_beta(TRAYON r); |
101 |
+ |
static double prob_beta_alpha(TRAYON r); |
102 |
+ |
static double prob_gamma_alpha(TRAYON r); |
103 |
+ |
static void v_par(FVECT v, FVECT v_out); |
104 |
+ |
static void v_per(FVECT v, FVECT v_out); |
105 |
+ |
static TRAYON transalphabeta(TRAYON r_initial); |
106 |
+ |
static TRAYON transbetaalpha(TRAYON r_initial); |
107 |
+ |
static TRAYON transalphagamma(TRAYON r_initial); |
108 |
+ |
static TRAYON transgammaalpha(TRAYON r_initial); |
109 |
+ |
static int compare(TRAYON r1, TRAYON r2, double marge); |
110 |
+ |
static void sortie(TRAYON r); |
111 |
+ |
static void trigo(TRAYON r); |
112 |
+ |
static TRAYON reflexion(TRAYON r_incident); |
113 |
+ |
static TRAYON transmission(TRAYON r_incident); |
114 |
+ |
static void trace_rayon(TRAYON r_incident); |
115 |
+ |
static void inverser(TRAYON *r1, TRAYON *r2); |
116 |
+ |
static void setprism(void); |
117 |
+ |
static double l_get_val(char *nm); |
118 |
|
|
119 |
+ |
|
120 |
|
/* Definition des routines */ |
121 |
|
|
122 |
|
#define term(a,b) a/Sqrt(a*a+b*b) |
123 |
< |
static |
124 |
< |
prepare_matrices() |
123 |
> |
static void |
124 |
> |
prepare_matrices(void) |
125 |
|
{ |
126 |
< |
/* preparation des matrices de changement de bases */ |
126 |
> |
/* preparation des matrices de changement de bases */ |
127 |
|
|
128 |
< |
matb[0][0] = matbt[0][0] = matb[1][1] = matbt[1][1] = term(prism.a,prism.d); |
129 |
< |
matb[1][0] = matbt[0][1] = term(-prism.d,prism.a); |
130 |
< |
matb[0][1] = matbt[1][0] = term(prism.d,prism.a); |
131 |
< |
matc[0][0] = matct[0][0] = matc[1][1] = matct[1][1] = term(prism.b,prism.d); |
132 |
< |
matc[1][0] = matct[0][1] = term(prism.d,prism.b); |
133 |
< |
matc[0][1] = matct[1][0] = term(-prism.d,prism.b); |
134 |
< |
return; |
128 |
> |
matb[0][0] = matbt[0][0] = matb[1][1] = matbt[1][1] = term(prism.a,prism.d); |
129 |
> |
matb[1][0] = matbt[0][1] = term(-prism.d,prism.a); |
130 |
> |
matb[0][1] = matbt[1][0] = term(prism.d,prism.a); |
131 |
> |
matc[0][0] = matct[0][0] = matc[1][1] = matct[1][1] = term(prism.b,prism.d); |
132 |
> |
matc[1][0] = matct[0][1] = term(prism.d,prism.b); |
133 |
> |
matc[0][1] = matct[1][0] = term(-prism.d,prism.b); |
134 |
> |
return; |
135 |
|
} |
136 |
|
#undef term |
137 |
|
|
138 |
|
|
139 |
< |
static |
140 |
< |
tfm(mat,v_old,v_new) |
141 |
< |
MAT4 mat; |
142 |
< |
FVECT v_old,v_new; |
139 |
> |
static void |
140 |
> |
tfm( |
141 |
> |
MAT4 mat, |
142 |
> |
FVECT v_old, |
143 |
> |
FVECT v_new |
144 |
> |
) |
145 |
|
{ |
146 |
< |
/* passage d'un repere old au repere new par la matrice mat */ |
147 |
< |
FVECT v_temp; |
146 |
> |
/* passage d'un repere old au repere new par la matrice mat */ |
147 |
> |
FVECT v_temp; |
148 |
|
|
149 |
< |
multv3(v_temp,v_old,mat); |
150 |
< |
normalize(v_temp); |
151 |
< |
VCOPY(v_new,v_temp); |
152 |
< |
return; |
149 |
> |
multv3(v_temp,v_old,mat); |
150 |
> |
normalize(v_temp); |
151 |
> |
VCOPY(v_new,v_temp); |
152 |
> |
return; |
153 |
|
} |
154 |
|
|
155 |
|
#define A prism.a |
159 |
|
|
160 |
|
|
161 |
|
static double |
162 |
< |
prob_alpha_beta(r) |
163 |
< |
TRAYON r; |
162 |
> |
prob_alpha_beta( |
163 |
> |
TRAYON r |
164 |
> |
) |
165 |
|
{ |
166 |
< |
/* calcul de la probabilite de passage de alpha a beta */ |
167 |
< |
double prob,test; |
166 |
> |
/* calcul de la probabilite de passage de alpha a beta */ |
167 |
> |
double prob,test; |
168 |
|
|
169 |
< |
if ( X(r) != 0. ) |
169 |
> |
if ( X(r) != 0. ) |
170 |
|
{ |
171 |
< |
test = Y(r)/X(r); |
172 |
< |
if ( test > B/D ) prob = 1.; |
173 |
< |
else if ( test >= -A/D ) prob = (A+test*D)/(A+B); |
174 |
< |
else prob = 0.; |
171 |
> |
test = Y(r)/X(r); |
172 |
> |
if ( test > B/D ) prob = 1.; |
173 |
> |
else if ( test >= -A/D ) prob = (A+test*D)/(A+B); |
174 |
> |
else prob = 0.; |
175 |
|
} |
176 |
< |
else prob = 0.; |
177 |
< |
return prob; |
176 |
> |
else prob = 0.; |
177 |
> |
return prob; |
178 |
|
} |
179 |
|
|
180 |
|
|
181 |
|
static double |
182 |
< |
prob_beta_alpha(r) |
183 |
< |
TRAYON r; |
182 |
> |
prob_beta_alpha( |
183 |
> |
TRAYON r |
184 |
> |
) |
185 |
|
{ |
186 |
< |
/* calcul de la probabilite de passage de beta a aplha */ |
187 |
< |
double prob,test; |
186 |
> |
/* calcul de la probabilite de passage de beta a aplha */ |
187 |
> |
double prob,test; |
188 |
|
|
189 |
< |
if ( X(r) != 0. ) |
189 |
> |
if ( X(r) != 0. ) |
190 |
|
{ |
191 |
< |
test = Y(r)/X(r); |
192 |
< |
if ( test > B/D ) prob = (A+B)/(A+test*D); |
193 |
< |
else if ( test >= -A/D ) prob = 1.; |
194 |
< |
else prob = 0.; |
191 |
> |
test = Y(r)/X(r); |
192 |
> |
if ( test > B/D ) prob = (A+B)/(A+test*D); |
193 |
> |
else if ( test >= -A/D ) prob = 1.; |
194 |
> |
else prob = 0.; |
195 |
|
} |
196 |
< |
else prob = 0.; |
197 |
< |
return prob; |
196 |
> |
else prob = 0.; |
197 |
> |
return prob; |
198 |
|
} |
199 |
|
|
200 |
|
|
201 |
< |
double prob_gamma_alpha(r) |
202 |
< |
TRAYON r; |
201 |
> |
static double |
202 |
> |
prob_gamma_alpha( |
203 |
> |
TRAYON r |
204 |
> |
) |
205 |
|
{ |
206 |
< |
/* calcul de la probabilite de passage de gamma a alpha */ |
207 |
< |
double prob,test; |
206 |
> |
/* calcul de la probabilite de passage de gamma a alpha */ |
207 |
> |
double prob,test; |
208 |
|
|
209 |
< |
if ( X(r) != 0. ) |
209 |
> |
if ( X(r) != 0. ) |
210 |
|
{ |
211 |
< |
test = Y(r)/X(r); |
212 |
< |
if ( test > B/D ) prob = 0.; |
213 |
< |
else if ( test >= -A/D ) prob = 1.; |
214 |
< |
else prob = (A+B)/(B-test*D); |
211 |
> |
test = Y(r)/X(r); |
212 |
> |
if ( test > B/D ) prob = 0.; |
213 |
> |
else if ( test >= -A/D ) prob = 1.; |
214 |
> |
else prob = (A+B)/(B-test*D); |
215 |
|
} |
216 |
< |
else prob = 0.; |
217 |
< |
return prob; |
216 |
> |
else prob = 0.; |
217 |
> |
return prob; |
218 |
|
} |
219 |
|
|
220 |
|
#undef A |
223 |
|
#undef D |
224 |
|
|
225 |
|
|
226 |
< |
static |
227 |
< |
v_par(v,v_out) |
228 |
< |
FVECT v,v_out; |
226 |
> |
static void |
227 |
> |
v_par( |
228 |
> |
FVECT v, |
229 |
> |
FVECT v_out |
230 |
> |
) |
231 |
|
/* calcule le vecteur par au plan d'incidence lie a v */ |
232 |
|
{ |
233 |
< |
FVECT v_temp; |
234 |
< |
double det; |
233 |
> |
FVECT v_temp; |
234 |
> |
double det; |
235 |
|
|
236 |
< |
det = Sqrt( (YY(v)*YY(v)+ZZ(v)*ZZ(v))*(YY(v)*YY(v)+ZZ(v)*ZZ(v))+ |
237 |
< |
(XX(v)*XX(v)*YY(v)*YY(v))+(XX(v)*XX(v)*ZZ(v)*ZZ(v)) ); |
238 |
< |
XX(v_temp) = (YY(v)*YY(v)+ZZ(v)*ZZ(v))/det; |
239 |
< |
YY(v_temp) = -( XX(v)*YY(v) )/det; |
240 |
< |
ZZ(v_temp) = -( XX(v)*ZZ(v) )/det; |
241 |
< |
VCOPY(v_out,v_temp); |
242 |
< |
return; |
236 |
> |
det = Sqrt( (YY(v)*YY(v)+ZZ(v)*ZZ(v))*(YY(v)*YY(v)+ZZ(v)*ZZ(v))+ |
237 |
> |
(XX(v)*XX(v)*YY(v)*YY(v))+(XX(v)*XX(v)*ZZ(v)*ZZ(v)) ); |
238 |
> |
XX(v_temp) = (YY(v)*YY(v)+ZZ(v)*ZZ(v))/det; |
239 |
> |
YY(v_temp) = -( XX(v)*YY(v) )/det; |
240 |
> |
ZZ(v_temp) = -( XX(v)*ZZ(v) )/det; |
241 |
> |
VCOPY(v_out,v_temp); |
242 |
> |
return; |
243 |
|
} |
244 |
|
|
245 |
|
|
246 |
< |
static |
247 |
< |
v_per(v,v_out) |
248 |
< |
FVECT v,v_out; |
246 |
> |
static void |
247 |
> |
v_per( |
248 |
> |
FVECT v, |
249 |
> |
FVECT v_out |
250 |
> |
) |
251 |
|
/* calcule le vecteur perp au plan d'incidence lie a v */ |
252 |
|
{ |
253 |
< |
FVECT v_temp; |
254 |
< |
double det; |
253 |
> |
FVECT v_temp; |
254 |
> |
double det; |
255 |
|
|
256 |
< |
det = Sqrt( (ZZ(v)*ZZ(v)+YY(v)*YY(v)) ); |
257 |
< |
XX(v_temp) = 0.; |
258 |
< |
YY(v_temp) = -ZZ(v)/det; |
259 |
< |
ZZ(v_temp) = YY(v)/det; |
260 |
< |
VCOPY(v_out,v_temp); |
261 |
< |
return; |
256 |
> |
det = Sqrt( (ZZ(v)*ZZ(v)+YY(v)*YY(v)) ); |
257 |
> |
XX(v_temp) = 0.; |
258 |
> |
YY(v_temp) = -ZZ(v)/det; |
259 |
> |
ZZ(v_temp) = YY(v)/det; |
260 |
> |
VCOPY(v_out,v_temp); |
261 |
> |
return; |
262 |
|
} |
263 |
|
|
264 |
|
|
265 |
|
static TRAYON |
266 |
< |
transalphabeta(r_initial) |
266 |
> |
transalphabeta( |
267 |
> |
TRAYON r_initial |
268 |
> |
) |
269 |
|
/* transforme le rayon r_initial de la base associee a alpha dans |
270 |
|
la base associee a beta */ |
239 |
– |
TRAYON r_initial; |
271 |
|
{ |
272 |
< |
TRAYON r_final; |
273 |
< |
FVECT vpar_temp1,vpar_temp2,vper_temp1,vper_temp2; |
272 |
> |
TRAYON r_final; |
273 |
> |
FVECT vpar_temp1,vpar_temp2,vper_temp1,vper_temp2; |
274 |
|
|
275 |
< |
r_final = r_initial; |
276 |
< |
alpha_beta(r_initial.v,r_final.v); |
277 |
< |
if ((Y(r_initial) != 0. || Z(r_initial) != 0.)&&(Y(r_final) !=0. || Z(r_final)!= 0.)) |
278 |
< |
{ |
279 |
< |
v_par(r_initial.v,vpar_temp1); |
280 |
< |
alpha_beta(vpar_temp1,vpar_temp1); |
281 |
< |
v_per(r_initial.v,vper_temp1); |
282 |
< |
alpha_beta(vper_temp1,vper_temp1); |
283 |
< |
v_par(r_final.v,vpar_temp2); |
284 |
< |
v_per(r_final.v,vper_temp2); |
285 |
< |
r_final.ppar1 = (r_initial.ppar1*fdot(vpar_temp1,vpar_temp2))+ |
286 |
< |
(r_initial.pper1*fdot(vper_temp1,vpar_temp2)); |
287 |
< |
r_final.pper1 = (r_initial.ppar1*fdot(vpar_temp1,vper_temp2))+ |
288 |
< |
(r_initial.pper1*fdot(vper_temp1,vper_temp2)); |
289 |
< |
r_final.ppar2 = (r_initial.ppar2*fdot(vpar_temp1,vpar_temp2))+ |
290 |
< |
(r_initial.pper2*fdot(vper_temp1,vpar_temp2)); |
291 |
< |
r_final.pper2 = (r_initial.ppar2*fdot(vpar_temp1,vper_temp2))+ |
292 |
< |
(r_initial.pper2*fdot(vper_temp1,vper_temp2)); |
293 |
< |
} |
294 |
< |
return r_final; |
275 |
> |
r_final = r_initial; |
276 |
> |
alpha_beta(r_initial.v,r_final.v); |
277 |
> |
if ((Y(r_initial) != 0. || Z(r_initial) != 0.)&&(Y(r_final) !=0. || Z(r_final)!= 0.)) |
278 |
> |
{ |
279 |
> |
v_par(r_initial.v,vpar_temp1); |
280 |
> |
alpha_beta(vpar_temp1,vpar_temp1); |
281 |
> |
v_per(r_initial.v,vper_temp1); |
282 |
> |
alpha_beta(vper_temp1,vper_temp1); |
283 |
> |
v_par(r_final.v,vpar_temp2); |
284 |
> |
v_per(r_final.v,vper_temp2); |
285 |
> |
r_final.ppar1 = (r_initial.ppar1*fdot(vpar_temp1,vpar_temp2))+ |
286 |
> |
(r_initial.pper1*fdot(vper_temp1,vpar_temp2)); |
287 |
> |
r_final.pper1 = (r_initial.ppar1*fdot(vpar_temp1,vper_temp2))+ |
288 |
> |
(r_initial.pper1*fdot(vper_temp1,vper_temp2)); |
289 |
> |
r_final.ppar2 = (r_initial.ppar2*fdot(vpar_temp1,vpar_temp2))+ |
290 |
> |
(r_initial.pper2*fdot(vper_temp1,vpar_temp2)); |
291 |
> |
r_final.pper2 = (r_initial.ppar2*fdot(vpar_temp1,vper_temp2))+ |
292 |
> |
(r_initial.pper2*fdot(vper_temp1,vper_temp2)); |
293 |
> |
} |
294 |
> |
return r_final; |
295 |
|
} |
296 |
|
|
297 |
|
|
298 |
|
static TRAYON |
299 |
< |
transbetaalpha(r_initial) |
300 |
< |
/* transforme le rayon r_initial de la base associee a beta dans |
301 |
< |
la base associee a alpha */ |
271 |
< |
TRAYON r_initial; |
299 |
> |
transbetaalpha( |
300 |
> |
TRAYON r_initial |
301 |
> |
) |
302 |
|
{ |
303 |
< |
TRAYON r_final; |
304 |
< |
FVECT vpar_temp1,vpar_temp2,vper_temp1,vper_temp2; |
303 |
> |
/* transforme le rayon r_initial de la base associee a beta dans |
304 |
> |
la base associee a alpha */ |
305 |
> |
TRAYON r_final; |
306 |
> |
FVECT vpar_temp1,vpar_temp2,vper_temp1,vper_temp2; |
307 |
|
|
308 |
< |
r_final = r_initial; |
309 |
< |
beta_alpha(r_initial.v,r_final.v); |
310 |
< |
if ((Y(r_initial) != 0. || Z(r_initial) != 0. )&&(Y(r_final) != 0. || Z(r_final)!= 0.)) |
311 |
< |
{ |
312 |
< |
v_par(r_initial.v,vpar_temp1); |
313 |
< |
beta_alpha(vpar_temp1,vpar_temp1); |
314 |
< |
v_per(r_initial.v,vper_temp1); |
315 |
< |
beta_alpha(vper_temp1,vper_temp1); |
316 |
< |
v_par(r_final.v,vpar_temp2); |
317 |
< |
v_per(r_final.v,vper_temp2); |
318 |
< |
r_final.ppar1 = (r_initial.ppar1*fdot(vpar_temp1,vpar_temp2))+ |
319 |
< |
(r_initial.pper1*fdot(vper_temp1,vpar_temp2)); |
320 |
< |
r_final.pper1 = (r_initial.ppar1*fdot(vpar_temp1,vper_temp2))+ |
321 |
< |
(r_initial.pper1*fdot(vper_temp1,vper_temp2)); |
322 |
< |
r_final.ppar2 = (r_initial.ppar2*fdot(vpar_temp1,vpar_temp2))+ |
323 |
< |
(r_initial.pper2*fdot(vper_temp1,vpar_temp2)); |
324 |
< |
r_final.pper2 = (r_initial.ppar2*fdot(vpar_temp1,vper_temp2))+ |
325 |
< |
(r_initial.pper2*fdot(vper_temp1,vper_temp2)); |
308 |
> |
r_final = r_initial; |
309 |
> |
beta_alpha(r_initial.v,r_final.v); |
310 |
> |
if ((Y(r_initial) != 0. || Z(r_initial) != 0. )&&(Y(r_final) != 0. || Z(r_final)!= 0.)) |
311 |
> |
{ |
312 |
> |
v_par(r_initial.v,vpar_temp1); |
313 |
> |
beta_alpha(vpar_temp1,vpar_temp1); |
314 |
> |
v_per(r_initial.v,vper_temp1); |
315 |
> |
beta_alpha(vper_temp1,vper_temp1); |
316 |
> |
v_par(r_final.v,vpar_temp2); |
317 |
> |
v_per(r_final.v,vper_temp2); |
318 |
> |
r_final.ppar1 = (r_initial.ppar1*fdot(vpar_temp1,vpar_temp2))+ |
319 |
> |
(r_initial.pper1*fdot(vper_temp1,vpar_temp2)); |
320 |
> |
r_final.pper1 = (r_initial.ppar1*fdot(vpar_temp1,vper_temp2))+ |
321 |
> |
(r_initial.pper1*fdot(vper_temp1,vper_temp2)); |
322 |
> |
r_final.ppar2 = (r_initial.ppar2*fdot(vpar_temp1,vpar_temp2))+ |
323 |
> |
(r_initial.pper2*fdot(vper_temp1,vpar_temp2)); |
324 |
> |
r_final.pper2 = (r_initial.ppar2*fdot(vpar_temp1,vper_temp2))+ |
325 |
> |
(r_initial.pper2*fdot(vper_temp1,vper_temp2)); |
326 |
|
|
327 |
< |
} |
328 |
< |
return r_final; |
327 |
> |
} |
328 |
> |
return r_final; |
329 |
|
} |
330 |
|
|
331 |
|
|
332 |
|
static TRAYON |
333 |
< |
transalphagamma(r_initial) |
333 |
> |
transalphagamma( |
334 |
> |
TRAYON r_initial |
335 |
> |
) |
336 |
|
/* transforme le rayon r_initial de la base associee a alpha dans |
337 |
|
la base associee a gamma */ |
304 |
– |
TRAYON r_initial; |
338 |
|
{ |
339 |
< |
TRAYON r_final; |
340 |
< |
FVECT vpar_temp1,vpar_temp2,vper_temp1,vper_temp2; |
339 |
> |
TRAYON r_final; |
340 |
> |
FVECT vpar_temp1,vpar_temp2,vper_temp1,vper_temp2; |
341 |
|
|
342 |
< |
r_final = r_initial; |
343 |
< |
alpha_gamma(r_initial.v,r_final.v); |
344 |
< |
if (( Y(r_initial) != 0. || Z(r_initial) != 0. )&&(Y(r_final)!= 0. || Z(r_final) !=0.)) |
345 |
< |
{ |
346 |
< |
v_par(r_initial.v,vpar_temp1); |
347 |
< |
alpha_gamma(vpar_temp1,vpar_temp1); |
348 |
< |
v_per(r_initial.v,vper_temp1); |
349 |
< |
alpha_gamma(vper_temp1,vper_temp1); |
350 |
< |
v_par(r_final.v,vpar_temp2); |
351 |
< |
v_per(r_final.v,vper_temp2); |
352 |
< |
r_final.ppar1 = (r_initial.ppar1*fdot(vpar_temp1,vpar_temp2))+ |
353 |
< |
(r_initial.pper1*fdot(vper_temp1,vpar_temp2)); |
354 |
< |
r_final.pper1 = (r_initial.ppar1*fdot(vpar_temp1,vper_temp2))+ |
355 |
< |
(r_initial.pper1*fdot(vper_temp1,vper_temp2)); |
356 |
< |
r_final.ppar2 = (r_initial.ppar2*fdot(vpar_temp1,vpar_temp2))+ |
357 |
< |
(r_initial.pper2*fdot(vper_temp1,vpar_temp2)); |
358 |
< |
r_final.pper2 = (r_initial.ppar2*fdot(vpar_temp1,vper_temp2))+ |
359 |
< |
(r_initial.pper2*fdot(vper_temp1,vper_temp2)); |
342 |
> |
r_final = r_initial; |
343 |
> |
alpha_gamma(r_initial.v,r_final.v); |
344 |
> |
if (( Y(r_initial) != 0. || Z(r_initial) != 0. )&&(Y(r_final)!= 0. || Z(r_final) !=0.)) |
345 |
> |
{ |
346 |
> |
v_par(r_initial.v,vpar_temp1); |
347 |
> |
alpha_gamma(vpar_temp1,vpar_temp1); |
348 |
> |
v_per(r_initial.v,vper_temp1); |
349 |
> |
alpha_gamma(vper_temp1,vper_temp1); |
350 |
> |
v_par(r_final.v,vpar_temp2); |
351 |
> |
v_per(r_final.v,vper_temp2); |
352 |
> |
r_final.ppar1 = (r_initial.ppar1*fdot(vpar_temp1,vpar_temp2))+ |
353 |
> |
(r_initial.pper1*fdot(vper_temp1,vpar_temp2)); |
354 |
> |
r_final.pper1 = (r_initial.ppar1*fdot(vpar_temp1,vper_temp2))+ |
355 |
> |
(r_initial.pper1*fdot(vper_temp1,vper_temp2)); |
356 |
> |
r_final.ppar2 = (r_initial.ppar2*fdot(vpar_temp1,vpar_temp2))+ |
357 |
> |
(r_initial.pper2*fdot(vper_temp1,vpar_temp2)); |
358 |
> |
r_final.pper2 = (r_initial.ppar2*fdot(vpar_temp1,vper_temp2))+ |
359 |
> |
(r_initial.pper2*fdot(vper_temp1,vper_temp2)); |
360 |
|
|
361 |
< |
} |
362 |
< |
return r_final; |
361 |
> |
} |
362 |
> |
return r_final; |
363 |
|
} |
364 |
|
|
365 |
|
|
366 |
|
static TRAYON |
367 |
< |
transgammaalpha(r_initial) |
367 |
> |
transgammaalpha( |
368 |
> |
TRAYON r_initial |
369 |
> |
) |
370 |
|
/* transforme le rayon r_initial de la base associee a gamma dans |
371 |
|
la base associee a alpha */ |
337 |
– |
TRAYON r_initial; |
372 |
|
{ |
373 |
< |
TRAYON r_final; |
374 |
< |
FVECT vpar_temp1,vpar_temp2,vper_temp1,vper_temp2; |
373 |
> |
TRAYON r_final; |
374 |
> |
FVECT vpar_temp1,vpar_temp2,vper_temp1,vper_temp2; |
375 |
|
|
376 |
< |
r_final = r_initial; |
377 |
< |
gamma_alpha(r_initial.v,r_final.v); |
378 |
< |
if (( Y(r_initial) != 0. || Z(r_initial) != 0. )&&(Y(r_final) !=0. || Z(r_final) != 0.)) |
379 |
< |
{ |
380 |
< |
v_par(r_initial.v,vpar_temp1); |
381 |
< |
gamma_alpha(vpar_temp1,vpar_temp1); |
382 |
< |
v_per(r_initial.v,vper_temp1); |
383 |
< |
gamma_alpha(vper_temp1,vper_temp1); |
384 |
< |
v_par(r_final.v,vpar_temp2); |
385 |
< |
v_per(r_final.v,vper_temp2); |
386 |
< |
r_final.ppar1 = (r_initial.ppar1*fdot(vpar_temp1,vpar_temp2))+ |
387 |
< |
(r_initial.pper1*fdot(vper_temp1,vpar_temp2)); |
388 |
< |
r_final.pper1 = (r_initial.ppar1*fdot(vpar_temp1,vper_temp2))+ |
389 |
< |
(r_initial.pper1*fdot(vper_temp1,vper_temp2)); |
390 |
< |
r_final.ppar2 = (r_initial.ppar2*fdot(vpar_temp1,vpar_temp2))+ |
391 |
< |
(r_initial.pper2*fdot(vper_temp1,vpar_temp2)); |
392 |
< |
r_final.pper2 = (r_initial.ppar2*fdot(vpar_temp1,vper_temp2))+ |
393 |
< |
(r_initial.pper2*fdot(vper_temp1,vper_temp2)); |
394 |
< |
} |
395 |
< |
return r_final; |
376 |
> |
r_final = r_initial; |
377 |
> |
gamma_alpha(r_initial.v,r_final.v); |
378 |
> |
if (( Y(r_initial) != 0. || Z(r_initial) != 0. )&&(Y(r_final) !=0. || Z(r_final) != 0.)) |
379 |
> |
{ |
380 |
> |
v_par(r_initial.v,vpar_temp1); |
381 |
> |
gamma_alpha(vpar_temp1,vpar_temp1); |
382 |
> |
v_per(r_initial.v,vper_temp1); |
383 |
> |
gamma_alpha(vper_temp1,vper_temp1); |
384 |
> |
v_par(r_final.v,vpar_temp2); |
385 |
> |
v_per(r_final.v,vper_temp2); |
386 |
> |
r_final.ppar1 = (r_initial.ppar1*fdot(vpar_temp1,vpar_temp2))+ |
387 |
> |
(r_initial.pper1*fdot(vper_temp1,vpar_temp2)); |
388 |
> |
r_final.pper1 = (r_initial.ppar1*fdot(vpar_temp1,vper_temp2))+ |
389 |
> |
(r_initial.pper1*fdot(vper_temp1,vper_temp2)); |
390 |
> |
r_final.ppar2 = (r_initial.ppar2*fdot(vpar_temp1,vpar_temp2))+ |
391 |
> |
(r_initial.pper2*fdot(vper_temp1,vpar_temp2)); |
392 |
> |
r_final.pper2 = (r_initial.ppar2*fdot(vpar_temp1,vper_temp2))+ |
393 |
> |
(r_initial.pper2*fdot(vper_temp1,vper_temp2)); |
394 |
> |
} |
395 |
> |
return r_final; |
396 |
|
} |
397 |
|
|
398 |
|
|
399 |
|
|
400 |
|
static int |
401 |
< |
compare(r1,r2,marge) |
402 |
< |
TRAYON r1, r2; |
403 |
< |
double marge; |
404 |
< |
|
401 |
> |
compare( |
402 |
> |
TRAYON r1, |
403 |
> |
TRAYON r2, |
404 |
> |
double marge |
405 |
> |
) |
406 |
|
{ |
407 |
< |
double arctg1, arctg2; |
407 |
> |
double arctg1, arctg2; |
408 |
|
|
409 |
< |
arctg1 = atan2(Y(r1),X(r1)); |
410 |
< |
arctg2 = atan2(Y(r2),X(r2)); |
411 |
< |
if ((arctg1 - marge <= arctg2) && (arctg1 + marge >= arctg2)) return 1; |
412 |
< |
else return 0; |
409 |
> |
arctg1 = atan2(Y(r1),X(r1)); |
410 |
> |
arctg2 = atan2(Y(r2),X(r2)); |
411 |
> |
if ((arctg1 - marge <= arctg2) && (arctg1 + marge >= arctg2)) return 1; |
412 |
> |
else return 0; |
413 |
|
} |
414 |
|
|
415 |
|
|
416 |
|
|
417 |
|
|
418 |
< |
static |
419 |
< |
sortie(r) |
420 |
< |
TRAYON r; |
418 |
> |
static void |
419 |
> |
sortie( |
420 |
> |
TRAYON r |
421 |
> |
) |
422 |
|
{ |
423 |
< |
int i = 0; |
424 |
< |
int egalite = 0; |
423 |
> |
int i = 0; |
424 |
> |
int egalite = 0; |
425 |
|
|
426 |
|
|
427 |
< |
if(r.e > seuil) |
428 |
< |
{ |
429 |
< |
while (i < nbrayons && egalite == 0) |
430 |
< |
{ |
431 |
< |
raytemp = &ray[i]; |
432 |
< |
egalite = compare(r,*raytemp,tolerance); |
433 |
< |
if (egalite) raytemp->e = raytemp->e + r.e; |
434 |
< |
else i = i + 1; |
435 |
< |
} |
436 |
< |
if (egalite == 0) |
437 |
< |
{ |
438 |
< |
if (nbrayons == 0) ray = (TRAYON *)calloc(1,sizeof(TRAYON)); |
439 |
< |
else ray = (TRAYON *)realloc((void *)ray, (nbrayons+1)*(sizeof(TRAYON))); |
440 |
< |
if (ray == NULL) |
441 |
< |
error(SYSTEM, "out of memory in sortie\n"); |
442 |
< |
raytemp = &ray[nbrayons]; |
443 |
< |
raytemp->v[0] = X(r); |
444 |
< |
raytemp->v[1] = Y(r); |
445 |
< |
raytemp->v[2] = Z(r); |
446 |
< |
raytemp->e = r.e; |
447 |
< |
nbrayons++; |
448 |
< |
} |
449 |
< |
} |
450 |
< |
return; |
427 |
> |
if(r.e > seuil) |
428 |
> |
{ |
429 |
> |
while (i < nbrayons && egalite == 0) |
430 |
> |
{ |
431 |
> |
raytemp = &ray[i]; |
432 |
> |
egalite = compare(r,*raytemp,tolerance); |
433 |
> |
if (egalite) raytemp->e = raytemp->e + r.e; |
434 |
> |
else i = i + 1; |
435 |
> |
} |
436 |
> |
if (egalite == 0) |
437 |
> |
{ |
438 |
> |
if (nbrayons == 0) ray = (TRAYON *)calloc(1,sizeof(TRAYON)); |
439 |
> |
else ray = (TRAYON *)realloc((void *)ray, (nbrayons+1)*(sizeof(TRAYON))); |
440 |
> |
if (ray == NULL) |
441 |
> |
error(SYSTEM, "out of memory in sortie\n"); |
442 |
> |
raytemp = &ray[nbrayons]; |
443 |
> |
raytemp->v[0] = X(r); |
444 |
> |
raytemp->v[1] = Y(r); |
445 |
> |
raytemp->v[2] = Z(r); |
446 |
> |
raytemp->e = r.e; |
447 |
> |
nbrayons++; |
448 |
> |
} |
449 |
> |
} |
450 |
> |
return; |
451 |
|
} |
452 |
|
|
453 |
|
|
454 |
< |
static |
455 |
< |
trigo(r) |
456 |
< |
TRAYON r; |
454 |
> |
static void |
455 |
> |
trigo( |
456 |
> |
TRAYON r |
457 |
> |
) |
458 |
|
/* calcule les grandeurs trigonometriques relatives au rayon incident |
459 |
|
et le rapport entre les indices du milieu refracteur et incident */ |
460 |
|
{ |
461 |
< |
double det; |
462 |
< |
|
463 |
< |
det = Sqrt(X(r)*X(r)+Y(r)*Y(r)+Z(r)*Z(r)); |
464 |
< |
sinus = Sqrt(Y(r)*Y(r)+Z(r)*Z(r))/det; |
465 |
< |
cosinus = Sqrt(X(r)*X(r))/det; |
466 |
< |
if (r.n == 1.) rapport = prism.np * prism.np; |
467 |
< |
else rapport = 1./(prism.np * prism.np); |
468 |
< |
return; |
461 |
> |
double det; |
462 |
> |
|
463 |
> |
det = Sqrt(X(r)*X(r)+Y(r)*Y(r)+Z(r)*Z(r)); |
464 |
> |
sinus = Sqrt(Y(r)*Y(r)+Z(r)*Z(r))/det; |
465 |
> |
cosinus = Sqrt(X(r)*X(r))/det; |
466 |
> |
if (r.n == 1.) rapport = prism.np * prism.np; |
467 |
> |
else rapport = 1./(prism.np * prism.np); |
468 |
> |
return; |
469 |
|
} |
470 |
|
|
471 |
|
|
472 |
|
static TRAYON |
473 |
< |
reflexion(r_incident) |
474 |
< |
TRAYON r_incident; |
473 |
> |
reflexion( |
474 |
> |
TRAYON r_incident |
475 |
> |
) |
476 |
|
{ |
477 |
|
/* calcul du rayon reflechi par une face */ |
478 |
|
TRAYON r_reflechi; |
513 |
|
|
514 |
|
|
515 |
|
static TRAYON |
516 |
< |
transmission(r_incident) |
517 |
< |
TRAYON r_incident; |
516 |
> |
transmission( |
517 |
> |
TRAYON r_incident |
518 |
> |
) |
519 |
|
{ |
520 |
|
/* calcul du rayon refracte par une face */ |
521 |
|
TRAYON r_transmis; |
563 |
|
if ( r_suite.e > seuil ) trace_rayon(r_suite) |
564 |
|
|
565 |
|
|
566 |
< |
static |
567 |
< |
trace_rayon(r_incident) |
568 |
< |
TRAYON r_incident; |
566 |
> |
static void |
567 |
> |
trace_rayon( |
568 |
> |
TRAYON r_incident |
569 |
> |
) |
570 |
|
{ |
571 |
|
/* trace le rayon donne */ |
572 |
|
TRAYON r_reflechi,r_transmis,r_suite; |
656 |
|
|
657 |
|
#undef ensuite |
658 |
|
|
659 |
< |
static |
660 |
< |
inverser(r1,r2) |
661 |
< |
TRAYON *r1,*r2; |
662 |
< |
|
659 |
> |
static void |
660 |
> |
inverser( |
661 |
> |
TRAYON *r1, |
662 |
> |
TRAYON *r2 |
663 |
> |
) |
664 |
|
{ |
665 |
< |
TRAYON temp; |
666 |
< |
temp = *r1; |
667 |
< |
*r1 = *r2; |
668 |
< |
*r2 = temp; |
665 |
> |
TRAYON temp; |
666 |
> |
temp = *r1; |
667 |
> |
*r1 = *r2; |
668 |
> |
*r2 = temp; |
669 |
|
} |
670 |
|
|
671 |
|
|
672 |
|
|
673 |
< |
static |
674 |
< |
setprism() |
673 |
> |
static void |
674 |
> |
setprism(void) |
675 |
|
{ |
676 |
|
double d; |
677 |
|
TRAYON r_initial,rsource; |
678 |
< |
int i,j,k; |
678 |
> |
int i,j; |
679 |
|
|
680 |
|
prismclock = eclock; |
681 |
|
r_initial.ppar1 = r_initial.pper2 = 1.; |
769 |
|
|
770 |
|
|
771 |
|
static double |
772 |
< |
l_get_val() |
773 |
< |
|
772 |
> |
l_get_val( |
773 |
> |
char *nm |
774 |
> |
) |
775 |
|
{ |
776 |
|
int val, dir, i, trouve, curseur; |
777 |
|
int nb; |
778 |
|
double valeur; |
779 |
< |
TRAYON *rayt, raynull; |
779 |
> |
TRAYON *rayt=NULL, raynull; |
780 |
|
|
781 |
|
if (prismclock < 0 || prismclock < eclock) setprism(); |
782 |
|
if (bidon == BADVAL) { |
823 |
|
} |
824 |
|
|
825 |
|
|
826 |
< |
setprismfuncs() |
826 |
> |
void |
827 |
> |
setprismfuncs(void) /* declared in func.h */ |
828 |
|
{ |
829 |
|
funset("fprism_val", 3, '=', l_get_val); |
830 |
|
} |