97a9fa22a38ca36265b783f2027c1e003b3dcc14
[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 <TFile.h>
26 #include <TSystem.h>
27 #include <TVector.h>
28
29 #include "AliFieldMap.h"
30 #include "AliMagFMaps.h"
31
32
33 ClassImp(AliMagFMaps)
34
35 //_______________________________________________________________________
36 AliMagFMaps::AliMagFMaps():
37   fSolenoid(0),
38   fSolenoidUser(0.),
39   fL3Option(0),
40   fFieldRead(0)
41 {
42   //
43   // Default constructor
44   //
45   fFieldMap[0] = fFieldMap[1] = fFieldMap[2] = 0;
46 }
47
48 //_______________________________________________________________________
49 AliMagFMaps::AliMagFMaps(const char *name, const char *title, const Int_t integ, 
50                          const Float_t factor, const Float_t fmax, const Int_t map, 
51                          const Int_t l3):
52   AliMagF(name,title,integ,factor,fmax),
53   fSolenoid(0),
54   fSolenoidUser(0),
55   fL3Option(l3),
56   fFieldRead(0)
57 {
58   //
59   // Standard constructor
60   //
61   fType         = kConMesh;
62   fFieldMap[0]  = 0;
63   fMap          = map;
64   fL3Option     = l3;
65   ReadField();
66   fFieldRead = 1;
67   //
68   // Don't replicate field information in gAlice
69   for (Int_t i = 0; i < 3; i++)  fFieldMap[i]->SetWriteEnable(0);
70   //
71 }
72
73 //_______________________________________________________________________
74 AliMagFMaps::AliMagFMaps(const AliMagFMaps &magf):
75   AliMagF(magf),
76   fSolenoid(0),
77   fL3Option(0),
78   fFieldRead(0)
79 {
80   //
81   // Copy constructor
82   //
83   magf.Copy(*this);
84 }
85
86 //_______________________________________________________________________
87 AliMagFMaps::~AliMagFMaps()
88 {
89   //
90   //  Destructor
91   //
92   delete fFieldMap[0];
93   delete fFieldMap[1];
94   delete fFieldMap[2];    
95 }
96
97 //_______________________________________________________________________
98 void AliMagFMaps::ReadField()
99 {
100   //  Read Field Map from file
101   //
102   //  don't read twice
103   //
104   if (fFieldRead) return;
105   fFieldRead = 1;
106   //    
107   char* fname;
108   TFile* file = 0;
109   if (fMap == k2kG) {
110       if (fL3Option) {
111           fSolenoid = 2.;
112           fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/L3B02.root");
113           file = new TFile(fname);
114           fFieldMap[0] = dynamic_cast<AliFieldMap*>(file->Get("L3B02"));
115           file->Close();
116           delete file;
117       }
118       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/DipB02.root");
119       file = new TFile(fname);
120       fFieldMap[1] = dynamic_cast<AliFieldMap*>(file->Get("DipB02"));
121       file->Close();
122       delete file;;
123       
124       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/ExtB02.root");
125       file = new TFile(fname);
126       fFieldMap[2] = dynamic_cast<AliFieldMap*>(file->Get("ExtB02"));
127       file->Close();
128       delete file;
129   } else if (fMap == k4kG) {
130       if (fL3Option) {
131           fSolenoid = 4.;
132           fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/L3B04.root");
133           file = new TFile(fname);
134           fFieldMap[0] = dynamic_cast<AliFieldMap*>(file->Get("L3B04"));
135           file->Close();
136           delete file;
137       }
138       
139       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/DipB04.root");
140       file = new TFile(fname);
141       fFieldMap[1] = dynamic_cast<AliFieldMap*>(file->Get("DipB04"));
142       file->Close();
143       delete file;
144       
145       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/ExtB04.root");
146       file = new TFile(fname);
147       fFieldMap[2] = dynamic_cast<AliFieldMap*>(file->Get("ExtB04"));
148       file->Close();
149       delete file;
150   } else if (fMap == k5kG) {
151       if (fL3Option) {
152           fSolenoid = 5.;
153           fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/L3B05.root");
154           file = new TFile(fname);
155           fFieldMap[0] = dynamic_cast<AliFieldMap*>(file->Get("L3B05"));
156           file->Close();
157           delete file;
158       }
159       
160       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/DipB05.root");
161       file = new TFile(fname);
162       fFieldMap[1] = dynamic_cast<AliFieldMap*>(file->Get("DipB05"));
163       file->Close();
164       delete file;
165       
166       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/ExtB05.root");
167       file = new TFile(fname);
168       fFieldMap[2] = dynamic_cast<AliFieldMap*>(file->Get("ExtB05"));
169       file->Close();
170       delete file;
171       
172   }
173   
174   if (!fL3Option) {
175       //
176       // Dummy L3 map
177       fFieldMap[0] = new AliFieldMap();
178       fFieldMap[0] -> SetLimits(-800., 800., -800., 800., -700., 700.);
179       switch(fMap) {
180       case k2kG:
181           fSolenoid = 2.;
182           break;
183       case k4kG:
184           fSolenoid = 4.;
185           break;
186       case k5kG:
187           fSolenoid = 5.;
188           break;
189       }
190       fSolenoidUser = fSolenoid;
191   }
192 }
193
194 //_______________________________________________________________________
195 Float_t AliMagFMaps::SolenoidField() const
196 {
197   //
198   // Returns max. L3 (solenoid) field strength 
199   // according to field map setting 
200   //
201   return fSolenoid;
202 }
203
204 //_______________________________________________________________________
205 void AliMagFMaps::Field(Float_t *x, Float_t *b)
206 {
207   //
208   // Method to calculate the magnetic field
209   //
210   // --- find the position in the grid ---
211   
212   if (!fFieldRead) ReadField();
213   
214   b[0]=b[1]=b[2]=0;
215   AliFieldMap* map = 0;
216   if (fFieldMap[0]->Inside(x[0], x[1], x[2])) {
217       map = fFieldMap[0];
218       if (!fL3Option) {
219       //
220       //     Constant L3 field, if this option was selected
221       //
222           b[2] = fSolenoidUser;
223           return;
224     }
225   } else if (fFieldMap[1]->Inside(x[0], x[1], x[2])) {
226     map = fFieldMap[1];
227   } else if (fFieldMap[2]->Inside(x[0], x[1], x[2])) {
228     map = fFieldMap[2];
229   }
230   
231   if(map){
232     map->Field(x,b);
233   } else {
234     //This is the ZDC part
235     Float_t rad2=x[0]*x[0]+x[1]*x[1];
236     if(x[2]>kCORBEG2 && x[2]<kCOREND2){
237         if(rad2<kCOR2RA2){
238             b[0] = kFCORN2;
239         }
240     } else if(x[2]>kZ1BEG && x[2]<kZ1END){  
241         if(rad2<kZ1RA2){
242             b[0] = -kG1*x[1];
243             b[1] = -kG1*x[0];
244         }
245     } else if(x[2]>kZ2BEG && x[2]<kZ2END){
246         if(rad2<kZ2RA2){
247             b[0] = kG1*x[1];
248             b[1] = kG1*x[0];
249         }
250     }
251     else if(x[2]>kZ3BEG && x[2]<kZ3END){
252         if(rad2<kZ3RA2){
253             b[0] = kG1*x[1];
254             b[1] = kG1*x[0];
255         }
256     }
257     else if(x[2]>kZ4BEG && x[2]<kZ4END){
258         if(rad2<kZ4RA2){
259             b[0] = -kG1*x[1];
260             b[1] = -kG1*x[0];
261         }
262     }
263     else if(x[2]>kD1BEG && x[2]<kD1END){ 
264         if(rad2<kD1RA2){
265             b[1] = -kFDIP;
266         }
267     }
268     else if(x[2]>kD2BEG && x[2]<kD2END){
269         if(((x[0]-kXCEN1D2)*(x[0]-kXCEN1D2)+(x[1]-kYCEN1D2)*(x[1]-kYCEN1D2))<kD2RA2
270            || ((x[0]-kXCEN2D2)*(x[0]-kXCEN2D2)+(x[1]-kYCEN2D2)*(x[1]-kYCEN2D2))<kD2RA2){
271             b[1] = kFDIP;
272         }
273     }
274   }
275   if(fFactor!=1) {
276     b[0]*=fFactor;
277     b[1]*=fFactor;
278     b[2]*=fFactor;
279   }
280 }
281
282 //_______________________________________________________________________
283 void AliMagFMaps::Copy(AliMagFMaps & /* magf */) const
284 {
285   //
286   // Copy *this onto magf -- Not implemented
287   //
288   Fatal("Copy","Not implemented!\n");
289 }
290
291 //_______________________________________________________________________
292 void AliMagFMaps::Streamer(TBuffer &R__b)
293 {
294   // Stream an object of class AliMagFMaps.
295   if (R__b.IsReading()) {
296     AliMagFMaps::Class()->ReadBuffer(R__b, this);
297     fFieldRead = 0;
298   } else {
299     AliMagFMaps::Class()->WriteBuffer(R__b, this);
300   }
301 }