7db65add35693b42ac9607d3d04a120f07433351
[u/mrichter/AliRoot.git] / STEER / AliMagFMaps.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 // Magnetic field composed by 3 maps: the L3 magnet, extended region, and
20 // dipole magnet
21 // Used in the configuration macros (macros/Config.C, etc.)
22 // Author: Andreas Morsch <andreas.morsch@cern.ch>
23 //------------------------------------------------------------------------
24
25 #include <TClass.h>
26 #include <TFile.h>
27 #include <TMath.h>
28 #include <TSystem.h>
29
30 #include "AliLog.h"
31 #include "AliFieldMap.h"
32 #include "AliMagFMaps.h"
33
34 ClassImp(AliMagFMaps)
35     
36
37 //_______________________________________________________________________
38 AliMagFMaps::AliMagFMaps():
39   fSolenoid(0),
40   fSolenoidUser(0.),
41   fL3Option(0)
42 {
43   //
44   // Default constructor
45   //
46   //
47   // Don't replicate field information in gAlice
48   fFieldMap[0] = fFieldMap[1] = fFieldMap[2] = 0;
49 }
50
51 //_______________________________________________________________________
52 AliMagFMaps::AliMagFMaps(const char *name, const char *title, Int_t integ, 
53                          Float_t factor, Float_t fmax, Int_t map, 
54                          Int_t l3):
55   AliMagFC(name,title,integ,factor,fmax),
56   fSolenoid(0),
57   fSolenoidUser(0),
58   fL3Option(l3)
59 {
60   //
61   // Standard constructor
62   //
63   fType         = kConMesh;
64   fFieldMap[0]  = 0;
65   fMap          = map;
66   fL3Option     = l3;
67   ReadField();
68   //
69   // Don't replicate field information in gAlice
70   for (Int_t i = 0; i < 3; i++)  fFieldMap[i]->SetWriteEnable(0);
71   //
72 }
73
74 //_______________________________________________________________________
75 AliMagFMaps::AliMagFMaps(const AliMagFMaps &magf):
76   AliMagFC(magf),
77   fSolenoid(0),
78   fSolenoidUser(0),
79   fL3Option(0)
80 {
81   //
82   // Copy constructor
83   //
84   magf.Copy(*this);
85 }
86
87 //_______________________________________________________________________
88 AliMagFMaps::~AliMagFMaps()
89 {
90   //
91   //  Destructor
92   //
93   delete fFieldMap[0];
94   delete fFieldMap[1];
95   delete fFieldMap[2];    
96 }
97
98 //_______________________________________________________________________
99 void AliMagFMaps::ReadField()
100 {
101   //  Read Field Map from file
102   //
103   //  don't read twice
104   //
105     if (!fReadField) return;
106     fReadField = kFALSE;
107   //    
108    
109
110   char* fname = 0;
111   TFile* file = 0;
112   if (fMap == k2kG) {
113       fSolenoid = 2.;
114       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/L3B02.root");
115       file = new TFile(fname);
116       fFieldMap[0] = dynamic_cast<AliFieldMap*>(file->Get("L3B02"));
117       file->Close();
118       delete file;
119       delete [] fname; // see documentation to expand filenam
120       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/DipB02.root");
121       file = new TFile(fname);
122       fFieldMap[1] = dynamic_cast<AliFieldMap*>(file->Get("DipB02"));
123       file->Close();
124       delete file;;
125       delete [] fname;
126       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/ExtB02.root");
127       file = new TFile(fname);
128       fFieldMap[2] = dynamic_cast<AliFieldMap*>(file->Get("ExtB02"));
129       file->Close();
130       delete file;
131       delete [] fname;
132   } else if (fMap == k4kG) {
133       fSolenoid = 4.;
134       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/L3B04.root");
135       file = new TFile(fname);
136       fFieldMap[0] = dynamic_cast<AliFieldMap*>(file->Get("L3B04"));
137       file->Close();
138       delete file;
139       delete [] fname;
140       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/DipB04.root");
141       file = new TFile(fname);
142       fFieldMap[1] = dynamic_cast<AliFieldMap*>(file->Get("DipB04"));
143       file->Close();
144       delete [] fname;
145       delete file;
146       
147       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/ExtB04.root");
148       file = new TFile(fname);
149       fFieldMap[2] = dynamic_cast<AliFieldMap*>(file->Get("ExtB04"));
150       file->Close();
151       delete [] fname;
152       delete file;
153   } else if (fMap == k5kG) {
154       fSolenoid = 5.;
155       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/L3B05.root");
156       file = new TFile(fname);
157       fFieldMap[0] = dynamic_cast<AliFieldMap*>(file->Get("L3B05"));
158       file->Close();
159       delete file;      
160       delete [] fname;
161       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/DipB05.root");
162       file = new TFile(fname);
163       fFieldMap[1] = dynamic_cast<AliFieldMap*>(file->Get("DipB05"));
164       file->Close();
165       delete file;
166       delete [] fname;
167       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/ExtB05.root");
168       file = new TFile(fname);
169       fFieldMap[2] = dynamic_cast<AliFieldMap*>(file->Get("ExtB05"));
170       file->Close();
171       delete file;
172       delete [] fname;
173   }
174
175   if (!fL3Option) fSolenoidUser = fSolenoid;
176
177 }
178
179 //_______________________________________________________________________
180 Float_t AliMagFMaps::SolenoidField() const
181 {
182   //
183   // Returns max. L3 (solenoid) field strength 
184   // according to field map setting 
185   //
186   return -Factor()*fSolenoid;
187 }
188
189 //_______________________________________________________________________
190 void AliMagFMaps::Field(const float *x, float *b) const
191 {
192   //
193   // Method to calculate the magnetic field
194   //
195   // --- find the position in the grid ---
196   
197
198     
199  //    if (!fFieldRead) ReadField();
200
201   //
202   // Field Maps have been calculated for the coordinate system in which 
203   // the Muon Spectrometer is placed at z > 0
204   // Transform coordinates corresponingly 
205   //
206   
207   b[0]=b[1]=b[2]=0;
208   Float_t xm[3];
209   xm[0] = - x[0];
210   xm[1] =   x[1];
211   xm[2] = - x[2];
212   
213   AliFieldMap* map = 0;
214   if (fFieldMap[0]->Inside(xm[0], xm[1], xm[2])) {
215       map = fFieldMap[0];
216       Float_t r = TMath::Sqrt(xm[0] * xm[0] + xm[1] * xm[1]);
217       
218       if (!fL3Option && TMath::Abs(xm[2]) < 370. && r < 550.) {
219       //
220       //     Constant L3 field , if this option was selected
221       //
222         b[2] = (- fSolenoid)*fFactor;
223           return;
224       } 
225   } else if (fFieldMap[1]->Inside(xm[0], xm[1], xm[2])) {
226     map = fFieldMap[1];
227   } else if (fFieldMap[2]->Inside(xm[0], xm[1], xm[2])) {
228     map = fFieldMap[2];
229   }
230   
231   if(map){
232       map->Field(xm,b);
233       b[0] = - b[0];
234       b[2] = - b[2];
235       
236       if(fFactor!=1) {
237           b[0]*=fFactor;
238           b[1]*=fFactor;
239           b[2]*=fFactor;
240       }
241   } else {
242       //This is the ZDC part
243       ZDCField(x, b);
244   }
245 }
246
247 //_______________________________________________________________________
248 void AliMagFMaps::Field(const double *x, double *b) const
249 {
250   //
251   // Method to calculate the magnetic field
252   //
253   // --- find the position in the grid ---
254   
255  //    if (!fFieldRead) ReadField();
256
257   //
258   // Field Maps have been calculated for the coordinate system in which 
259   // the Muon Spectrometer is placed at z > 0
260   // Transform coordinates corresponingly 
261   //
262   
263   b[0]=b[1]=b[2]=0;
264   double xm[3];
265   xm[0] = - x[0];
266   xm[1] =   x[1];
267   xm[2] = - x[2];
268   
269   AliFieldMap* map = 0;
270   if (fFieldMap[0]->Inside(xm[0], xm[1], xm[2])) {
271       map = fFieldMap[0];
272       Float_t r = TMath::Sqrt(xm[0] * xm[0] + xm[1] * xm[1]);
273       
274       if (!fL3Option && TMath::Abs(xm[2]) < 370. && r < 550.) {
275       //
276       //     Constant L3 field , if this option was selected
277       //
278         b[2] = (- fSolenoid)*fFactor;
279           return;
280       } 
281   } else if (fFieldMap[1]->Inside(xm[0], xm[1], xm[2])) {
282     map = fFieldMap[1];
283   } else if (fFieldMap[2]->Inside(xm[0], xm[1], xm[2])) {
284     map = fFieldMap[2];
285   }
286   
287   if(map){
288       map->Field(xm,b);
289       b[0] = - b[0];
290       b[2] = - b[2];
291       
292       if(fFactor!=1) {
293           b[0]*=fFactor;
294           b[1]*=fFactor;
295           b[2]*=fFactor;
296       }
297   } else {
298       //This is the ZDC part
299       ZDCField(x, b);
300   }
301 }
302
303 //_______________________________________________________________________
304 void AliMagFMaps::Copy(TObject & /* magf */) const
305 {
306   //
307   // Copy *this onto magf -- Not implemented
308   //
309   AliFatal("Not implemented!");
310 }
311
312 //_______________________________________________________________________
313 void AliMagFMaps::Streamer(TBuffer &R__b)
314 {
315   // Stream an object of class AliMagFMaps.
316   if (R__b.IsReading()) {
317    R__b.ReadClassBuffer(AliMagFMaps::Class(), this);
318    // 
319    fReadField = kTRUE;
320    ReadField();
321   } else {
322     R__b.WriteClassBuffer(AliMagFMaps::Class(), this);
323   }
324 }