Correction in OnionCone.
[u/mrichter/AliRoot.git] / STRUCT / AliALIFE.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /*
17 $Log$
18 Revision 1.2  2000/10/10 06:40:25  hristov
19 Bug fixes
20
21 Revision 1.1  2000/07/26 15:10:57  morsch
22 Helper class to write geometry in ALIFE format in parallel with Geant geometry definition.
23
24 */
25
26 #include <AliALIFE.h>
27
28 ClassImp(AliALIFE)
29
30     AliALIFE::AliALIFE(const char* file1, const char* file2)
31 {  
32 // Constructor
33     fNBodies = 0; 
34     fNVolumes= 0; 
35     fBodyFile   = file1;   // File for Fluka bodies
36     fVolumeFile = file2;
37     fFile1=fopen(fBodyFile,"w");
38     fFile2=fopen(fVolumeFile,"w");
39     BodyHeader();
40     VolumeHeader();
41     fDefaultVolume1 = "DEFAU1";
42     fDefaultVolume2 = "DEFAU2";    
43 }
44
45     AliALIFE::AliALIFE()
46 {
47 // Default constructor
48     fBodyFile   = "FlukaBody.inp";  
49     fVolumeFile = "FlukaVolume.inp";
50     fFile1=fopen(fBodyFile,"w");
51     fFile2=fopen(fVolumeFile,"w");
52     BodyHeader();
53     VolumeHeader();
54 }
55
56 void AliALIFE::BodyHeader()
57 {
58 // Write header for body definitions
59     fprintf(fFile1,"*\n");        
60     fprintf(fFile1,"GEOBEGIN                                                              COMBINAT\n");
61     fprintf(fFile1,"    0    0         AliRoot Generated\n");
62     fprintf(fFile1,"*\n");
63     fprintf(fFile1,"*== Body Definitions ================================================\n");
64     fprintf(fFile1,"*\n");
65                             
66 }
67
68
69 void AliALIFE::VolumeHeader()
70 {
71 // Write header for region (volume)  definitions
72     fprintf(fFile2,"REGIONS\n");
73     fprintf(fFile2,"*\n");
74     fprintf(fFile2,"*== Region Definitions =================================================\n");
75     fprintf(fFile2,"*\n");
76 }
77
78
79 void AliALIFE:: Cylinder(Float_t rmin, Float_t rmax, 
80                          Float_t zmin, Float_t zmax, 
81                          Float_t pos[3], 
82                          char* Material, char* Field, char* Cuts) 
83 {
84 // Simple cylinder
85 //
86 // Bodies
87 // ^^^^^^
88     char name1[5], name2[5], name3[5], name4[5];
89
90 //  outer radius
91     sprintf(name1, "R%4.4d", fNBodies++);
92     fprintf(fFile1,"%5s ZCC%10.3f%10.3f%10.3f\n",name1, pos[0], pos[1], rmax); // inner radius
93     sprintf(name2, "R%4.4d", fNBodies++);
94     fprintf(fFile1,"%5s ZCC%10.3f%10.3f%10.3f\n",name2, pos[0], pos[1], rmin); 
95 // z-max
96     sprintf(name3, "Z%4.4d", fNBodies++);
97     fprintf(fFile1,"%5s XYP%10.3f\n",name3, zmax); 
98 // z-min
99     sprintf(name4, "Z%4.4d", fNBodies++);
100     fprintf(fFile1,"%5s XYP%10.3f\n",name4, zmin); 
101 //
102 // Volumes
103 // ^^^^^^^
104
105     fprintf(fFile2,">%s:%s\n", Material, Field);
106     fprintf(fFile2,"EMFCUT=%s\n", Cuts);
107     fprintf(fFile2,"WW-FACTOR=%s\n", Cuts);
108     if (rmin >0) {
109         fprintf(fFile2,"+%s-%s+%s-%s\n", name1, name2, name3, name4);
110     } else {
111         fprintf(fFile2,"+%s+%s-%s\n", name1, name3, name4);
112     }
113     
114     fprintf(fFile2,"\n");    
115 }
116
117 void AliALIFE::OnionCylinder(Float_t* r, Int_t nr, Float_t zmin, Float_t zmax,
118                              Float_t pos[3], char** Materials, 
119                              char** Fields, char** Cuts) 
120 {
121 //
122 // Concentric cylinders
123 //
124
125 // Bodies
126 // ^^^^^^
127     char nameRin[6], nameRou[6], nameZin[6], nameZou[6];
128 // z-limits
129 // z-max
130     sprintf(nameZou, "Z%4.4d", fNBodies++);
131     nameZou[5]='\0';
132     
133     fprintf(fFile1,"%5s XYP%10.3f\n",nameZou, zmax); 
134 // z-min
135     sprintf(nameZin, "Z%4.4d", fNBodies++);
136     nameZin[5]='\0';
137     fprintf(fFile1,"%5s XYP%10.3f\n",nameZin, zmin); 
138
139 // inner radius
140     sprintf(nameRin, "R%4.4d", fNBodies++);
141     nameRin[5]='\0';
142     fprintf(fFile1,"%5s ZCC%10.3f%10.3f%10.3f\n",nameRin, pos[0], pos[1], r[0]); 
143
144     Int_t i;
145     for (i=1; i<nr; i++) {
146 //  outer radius
147         sprintf(nameRou, "R%4.4d", fNBodies++);
148         nameRou[5]='\0';
149         fprintf(fFile1,"%5s ZCC%10.3f%10.3f%10.3f\n",
150                 nameRou, pos[0], pos[1], r[i]);         
151 //
152 // Volumes
153 // ^^^^^^^
154         if (Fields && Cuts) {
155             fprintf(fFile2,">%s:%s\n", Materials[i-1], Fields[i-1]);
156             fprintf(fFile2,"EMFCUT=%s\n", Cuts[i-1]);
157             fprintf(fFile2,"WW-FACTOR=%s\n", Cuts[i-1]);
158         } else {
159             fprintf(fFile2,">%s:%s\n", Materials[i-1], "MF");
160             fprintf(fFile2,"EMFCUT=%s\n", "$UNSHIELDED");
161             fprintf(fFile2,"WW-FACTOR=%s\n", "$UNSHIELDED");
162         }
163         if (r[i-1] != 0.) {
164             fprintf(fFile2,"+%5s-%5s+%5s-%5s\n", nameZou, nameZin, nameRou, nameRin);
165         } else {
166             fprintf(fFile2,"+%5s-%5s+%5s\n", nameZou, nameZin, nameRou);
167         }
168         fprintf(fFile2,"\n");
169         strcpy(nameRin,nameRou);
170     }
171 }
172
173
174 void AliALIFE::Cone(Float_t rmin1, Float_t rmin2, 
175                      Float_t rmax1, Float_t rmax2,
176                      Float_t zmin, Float_t zmax, 
177                      Float_t pos[3], 
178                      char* Material, char* Field, char* Cuts) 
179 {
180 // Simple cone 
181
182 //
183 // Bodies
184 // ^^^^^^
185     char nameCou[6], nameCin[6];
186     Float_t d, r1, r2;
187     char nameZin[6], nameZou[6];
188 // z-max
189     sprintf(nameZou, "Z%4.4d", fNBodies++);
190     fprintf(fFile1,"%5s XYP%10.3f\n",nameZou, zmax); 
191 // z-min
192     sprintf(nameZin, "Z%4.4d", fNBodies++);
193     fprintf(fFile1,"%5s XYP%10.3f\n",nameZin, zmin); 
194     
195 //  outer radius
196     d=zmax-zmin;
197     if (rmax1 >= 0. && rmax2 >= 0.) {
198         if (rmax1!=rmax2) {
199             if (rmax1 > rmax2) {
200                 pos[2]=zmin;
201                 r1=rmax1;
202                 r2=rmax2;
203             } else {
204                 d=-d;
205                 pos[2]=zmax;
206                 r1=rmax2;
207                 r2=rmax1;
208             }
209             sprintf(nameCou, "C%4.4d", fNBodies++);
210             fprintf(fFile1,"%5s TRC%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f\n",
211                     nameCou, pos[0], pos[1], pos[2], 0., 0., d); 
212             fprintf(fFile1,"         %10.3f%10.3f\n",r1,r2);
213         } else {
214             sprintf(nameCou, "C%4.4d", fNBodies++);
215             fprintf(fFile1,"%5s ZCC%10.3f%10.3f%10.3f\n",
216                     nameCou, pos[0], pos[1], rmax1); 
217         } 
218     }else {
219         strcpy(nameCou,fDefaultVolume1);
220     }
221
222     
223     
224 // inner radius
225     if (rmin1 >= 0. && rmin2 >= 0.) {
226         if (rmin1!=rmin2) {
227             if (rmin1 != 0 && rmin2 !=0) {
228                 d=zmax-zmin;
229                 if (rmin1 > rmin2) {
230                     pos[2]=zmin;
231                     r1=rmin1;
232                     r2=rmin2;
233                 } else {
234                     pos[2]=zmax;
235                     r1=rmin2;
236                     r2=rmin1;
237                     d=-d;
238                 }
239                 sprintf(nameCin, "C%4.4d", fNBodies++);
240                 fprintf(fFile1,"%5s TRC%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f\n",
241                         nameCin, pos[0], pos[1], pos[2], 0., 0., d); 
242                 fprintf(fFile1,"         %10.3f%10.3f\n",r1,r2);
243             } 
244         } else {
245             sprintf(nameCin, "C%4.4d", fNBodies++);
246             fprintf(fFile1,"%5s ZCC%10.3f%10.3f%10.3f\n",
247                     nameCin, pos[0], pos[1], rmin1); 
248         } 
249     }else {
250         strcpy(nameCin,fDefaultVolume2);
251     }
252
253
254
255 //
256 // Volumes
257 // ^^^^^^^
258     fprintf(fFile2,">%s:%s\n", Material, Field);
259     fprintf(fFile2,"EMFCUT=%s\n", Cuts);
260     fprintf(fFile2,"WW-FACTOR=%s\n", Cuts);
261     if (rmin1 != 0 && rmin2 !=0) {
262         fprintf(fFile2,"+%s-%s+%s-%s\n", nameCou, nameCin, nameZou, nameZin);
263     } else {
264         fprintf(fFile2,"+%s+%s-%s\n", nameCou, nameZou, nameZin);
265     }
266     fprintf(fFile2,"\n");        
267 }
268
269 void AliALIFE::OnionCone (Float_t* r1, Float_t* r2, Int_t nr, 
270                           Float_t zmin, Float_t zmax,
271                           Float_t pos[3], char** Materials, char** Fields,
272                           char** Cuts) 
273 {
274 // Concentric cones
275 //
276     char nameCin[6], nameCou[6], nameZin[6], nameZou[6];
277 //
278 // Bodies
279 // ^^^^^^
280     Float_t ri1, ri2, ro1, ro2, d;
281
282 // z-max
283     sprintf(nameZou, "Z%4.4d", fNBodies++);
284     fprintf(fFile1,"%5s XYP%10.3f\n",nameZou, zmax); 
285 // z-min
286     sprintf(nameZin, "Z%4.4d", fNBodies++);
287     fprintf(fFile1,"%5s XYP%10.3f\n",nameZin, zmin); 
288      
289 //  inner radius
290     d=zmax-zmin;
291     Bool_t hasInner=kFALSE;
292     Bool_t isCylinder=kFALSE;
293     if (r1[0]>=0 && r2[0]>=0) {
294         if (r1[0] != 0. &&  r2[0] !=0.) {
295             if (r1[0]!=r2[0]) {
296                 hasInner=kTRUE;
297                 if (r1[0] > r2[0]) {
298                     pos[2]=zmin;
299                     ri1=r1[0];
300                     ri2=r2[0];
301                 } else {
302                     d=-d;
303                     pos[2]=zmax;
304                     ri1=r2[0];
305                     ri2=r1[0];
306                 }
307                 sprintf(nameCin, "C%4.4d", fNBodies++);
308                 nameCin[5]='\0';
309                 fprintf(fFile1,"%5s TRC%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f\n",
310                         nameCin, pos[0], pos[1], pos[2], 0., 0., d); 
311                 fprintf(fFile1,"         %10.3f%10.3f\n",ri1,ri2);
312             } else {
313                 isCylinder=kTRUE;
314                 sprintf(nameCin, "C%4.4d", fNBodies++);
315                 nameCin[5]='\0';
316                 fprintf(fFile1,"%5s ZCC%10.3f%10.3f%10.3f\n",
317                         nameCin, pos[0], pos[1], r1[0]); 
318             }
319         }
320     } else {
321         strcpy(nameCin,fDefaultVolume1);
322     }
323     
324     
325     
326 // outer radius
327     Int_t i;
328     for (i=1; i<nr; i++) {
329         if (r1[i] >= 0. && r2[i] >=0) {
330             if (r1[i]!=r2[i]) {
331                 d=zmax-zmin;
332                 if (r1[i] > r2[i]) {
333                     pos[2]=zmin;
334                     ro1=r1[i];
335                     ro2=r2[i];
336                 } else {
337                     pos[2]=zmax;
338                     ro1=r2[i];
339                     ro2=r1[i];
340                     d=-d;
341                 }
342                 sprintf(nameCou, "C%4.4d", fNBodies++);
343                 nameCou[5]='\0';
344                 fprintf(fFile1,"%5s TRC%10.3f%10.3f%10.3f%10.3f%10.3f%10.3f\n",
345                         nameCou, pos[0], pos[1], pos[2], 0., 0., d); 
346                 fprintf(fFile1,"         %10.3f%10.3f\n",ro1,ro2);
347             } else {
348                 isCylinder=kTRUE;
349                 sprintf(nameCou, "C%4.4d", fNBodies++);
350                 nameCou[5]='\0';
351                 fprintf(fFile1,"%5s ZCC%10.3f%10.3f%10.3f\n",
352                         nameCou, pos[0], pos[1], r1[i]); 
353             }
354         } else {
355             strcpy(nameCou,fDefaultVolume1);
356         }
357         
358 // Volumes
359 // ^^^^^^^
360         if (Fields && Cuts) {
361             fprintf(fFile2,">%s:%s\n", Materials[i-1], Fields[i-1]);
362             fprintf(fFile2,"EMFCUT=%s\n", Cuts[i-1]);
363             fprintf(fFile2,"WW-FACTOR=%s\n", Cuts[i-1]);
364         } else {
365             fprintf(fFile2,">%s:%s\n", Materials[i-1], "MF");
366             fprintf(fFile2,"EMFCUT=%s\n", "$UNSHIELDED");
367             fprintf(fFile2,"WW-FACTOR=%s\n", "$UNSHIELDED");
368         }
369         if (hasInner || isCylinder) {
370             if (!isCylinder) fprintf(fFile2,"+%5s-%5s+%5s-%5s\n", 
371                     nameCou, nameCin, nameZou, nameZin);
372             if (isCylinder)  fprintf(fFile2,"+%5s+%5s-%5s\n",
373                     nameCou, nameZou, nameZin); 
374         } else {
375             fprintf(fFile2,"+%5s\n", nameCou);
376             hasInner=kTRUE;
377         }
378         fprintf(fFile2,"\n");        
379         strcpy(nameCin,nameCou);
380     }
381 }
382
383
384 void AliALIFE::PolyCone(Float_t* rmin, Float_t* rmax, Float_t* z, 
385                         Int_t nz,
386                         Float_t pos[3], 
387                         char* Material, char* Field, char* Cuts) 
388 {
389 //
390 // Equivalent to the Geant3 PCON
391     Int_t i;
392     for (i=0; i<nz-1; i++) {
393 // skip (z_i = z_i+1)
394         if (z[i] == z[i+1]) continue;
395         Cone(rmin[i], rmin[i+1], rmax[i], rmax[i+1], z[i], z[i+1], pos,
396              Material, Field, Cuts);
397     }
398 }
399
400 void AliALIFE::OnionPolyCone(Float_t** r, Float_t* z,
401                              Int_t nr, Int_t nz,
402                              Float_t pos[3], 
403                              char** Materials, char** Fields, char** Cuts)
404 {
405 //
406 // Concentric PCONS
407     Int_t i, j;
408     for (i=0; i<nz-1; i++) {
409 // skip (z_i = z_i+1)
410         if (z[i] == z[i+1]) continue;
411         for (j=0; j<nr-1; j++) {
412             if (Fields && Cuts) {
413                 Cone(r[i][j], r[i+1][j], r[i][j+1], r[i+1][j+1], 
414                      z[i], z[i+1], pos,
415                      Materials[i], Fields[i], Cuts[i]);
416             } else {
417                 Cone(r[i][j], r[i+1][j], r[i][j+1], r[i+1][j+1], 
418                      z[i], z[i+1], pos,
419                      Materials[i]);
420             }
421         }
422     }
423 }
424
425 void AliALIFE::Comment(char* Comment)
426 {
427 // Insert comment line
428     fprintf(fFile1,"*%s\n", Comment);        
429     fprintf(fFile2,"*%s\n", Comment);        
430 }
431
432
433 void AliALIFE::Finish()
434 {
435 // Finish geometry definition
436     char s[BUFSIZ];
437     fprintf(fFile1,"      END\n");
438     fclose(fFile2);
439     fFile2=fopen(fVolumeFile,"r");
440     while (fgets(s, BUFSIZ, fFile2)) {
441         fputs(s,fFile1);
442     }
443     
444     fclose(fFile1);
445     fclose(fFile2);    
446 }
447
448
449
450
451
452
453
454
455
456
457