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