Transition to NewIO
[u/mrichter/AliRoot.git] / macros / plotField.C
1 void plotField(Int_t iField = 0)
2 {
3 //
4 //  iField = 0    2 kG  solenoid 
5 //           1    4 kG  solenoid
6 //           2    5 kG  solenoid
7 //  
8 //  load necessary libraries
9     gSystem->Load("$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET)/libminicern");
10     gSystem->Load("$(ROOTSYS)/lib/libPhysics");
11     gSystem->Load("$(ROOTSYS)/lib/libEG");
12     gSystem->Load("$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET)/libSTEER");
13     
14 //
15 //
16 //  create field map
17
18      AliMagFMaps* field = new AliMagFMaps(
19      "Maps","Maps", 
20      2, 1., 10., iField);
21
22 //     field-SetL3ConstField(1);
23      
24 //
25 //  get parameters
26      Float_t xMin, xMax, yMin, yMax, zMin, zMax;
27      Float_t dX, dY, dZ;
28      xMin = -350.;
29      xMax =  350.;
30      dX   = (field->FieldMap(0))->DelX();
31      yMin = -350.;
32      yMax =  350.;
33      dY   = (field->FieldMap(0))->DelY();
34      zMin =   250.;
35      zMax =  1450.;
36      dZ   = (field->FieldMap(0))->DelZ();
37      Int_t nx = (xMax-xMin)/dX;
38      Int_t ny = (yMax-yMin)/dY;
39      Int_t nz = (zMax-zMin)/dZ;
40      
41          
42 //
43 // create histogram
44      TH2F* hMap = new TH2F("hMap", "Field Map y-z", 
45                            nz, zMin, zMax, ny, yMin, yMax);
46      TH2F* hMap1 = new TH2F("hMap1", "Field Map y-z", 
47                            nz, zMin, zMax, ny, yMin, yMax);
48      TH1F* hZ   = new TH1F("hZ", "Field along Z", 
49                            nz, zMin, zMax);
50
51      TH2F* hVec = new TH2F("hVec", "Field Map y-z", 
52                            nz, zMin, zMax, ny, yMin, yMax);
53
54      TH2F* hCir = new TH2F("hCir", "Field Map y-z", 
55                            nz, zMin, zMax, ny, yMin, yMax);
56      Float_t bMax = 0.;
57      for (Int_t i = 0; i < nz; i++) {
58          for (Int_t j = 0; j < ny; j++) {
59              Float_t x[3];
60              Float_t b[3];
61              
62              x[2] = zMin + i * dZ;
63              x[1] = yMin + j * dY;
64              x[0] = 0.;
65              field->Field(x, b);
66              Float_t bb = TMath::Sqrt(b[0]*b[0]+b[1]*b[1]+b[2]*b[2]);
67              if (bb > bMax) bMax = bb;
68              hMap->Fill(x[2], x[1], bb);
69          }
70      }
71
72      for (Int_t i = 0; i < nz; i++) {
73          for (Int_t j = 0; j < ny; j++) {
74              x[2] = zMin + i * dZ;
75              x[1] = yMin + j * dY;
76              x[0] = 0.;
77              field->Field(x, b);
78              Float_t bb = TMath::Sqrt(b[0]*b[0]+b[1]*b[1]+b[2]*b[2]);
79              Float_t db = (bMax-bb)/bMax*100.;
80              
81              hMap1->Fill(x[2], x[1], db);
82          }
83      }
84      
85      for (Int_t i = 0; i < nz; i++) {
86          x[2] = zMin + i * dZ +dZ/2.;
87          x[1] = 0.;
88          x[0] = 0.;
89          field->Field(x, b);
90          Float_t bb = TMath::Sqrt(b[0]*b[0]+b[1]*b[1]+b[2]*b[2]);
91          hZ->Fill(x[2], bb);
92      }
93
94      TCanvas *c1 = new TCanvas("c1","Canvas 1",400,10,600,700);
95      hMap->Draw();
96      TCanvas *c2 = new TCanvas("c2","Canvas 2",400,10,600,700);
97      hZ->Draw();
98      TCanvas *c3 = new TCanvas("c3","Canvas 3",400,10,600,700);
99      hVec->Draw();
100      
101      Float_t scale1 = 0.9*TMath::Sqrt(dZ*dZ+dY*dY)/bMax/2.;
102      Float_t scale2 = 0.005/bMax;
103      TArrow* arrow;
104      
105      
106      for (Int_t i = 0; i < nz; i++) {
107          for (Int_t j = 0; j < nx; j++) {
108              x[2] = zMin + i * dZ + dZ/2.;
109              x[0] = xMin + j * dX + dX/2.;
110              x[1] = 0.;
111              field->Field(x, b);
112              Float_t bb = TMath::Sqrt(b[0]*b[0]+b[1]*b[1]+b[2]*b[2]);
113              b[2] *= scale1;
114              b[0] *= scale1;
115              Float_t width = 0.005+scale2*b[1];
116              arrow = new TArrow(x[2], x[0], x[2]+b[2], x[0]+b[0], width,"|>");
117              arrow->SetFillColor(1);
118              arrow->SetFillStyle(1001);
119              arrow->Draw();
120              c3->Modified();
121              c3->cd();
122          }
123      }
124
125      TCanvas *c4 = new TCanvas("c4","Canvas 4",400,10,600,700);
126      hCir->Draw();
127      for (Int_t i = 0; i < nz; i++) {
128          for (Int_t j = 0; j < ny; j++) {
129              x[2] = zMin + i * dZ + dZ/2.;
130              x[1] = yMin + j * dY + dY/2.;
131              x[0] = 0.;
132              field->Field(x, b);
133              Float_t bb = TMath::Sqrt(b[0]*b[0]+b[1]*b[1]+b[2]*b[2]);
134              TEllipse *ellipse; 
135              ellipse= new TEllipse(x[2], x[1], b[1]*scale1/4., b[1]*scale1/4.,0,360,0);
136              ellipse->Draw();
137              c4->Modified();
138              c4->cd();
139          }
140      }
141 }
142
143
144