Skip to content

Commit e706748

Browse files
njacazionjacazio
authored
[ALICE3] MagField: add short mag field, add visualization (#15700)
* Add short mag field macro
1 parent bff2b37 commit e706748

3 files changed

Lines changed: 183 additions & 29 deletions

File tree

‎Detectors/Upgrades/ALICE3/macros/ALICE3Field.C‎

Lines changed: 69 additions & 28 deletions
Original file line numberDiff line numberDiff line change
@@ -11,56 +11,97 @@
1111
//
1212
// Author: J. E. Munoz Mendez jesus.munoz@cern.ch
1313

14+
#include <functional>
15+
#include <cmath>
16+
#include <TCanvas.h>
17+
#include <TH2F.h>
18+
#include <TStyle.h>
19+
1420
std::function<void(const double*, double*)> field()
1521
{
1622
return [](const double* x, double* b) {
17-
double Rc;
18-
double R1;
19-
double R2;
20-
double B1;
21-
double B2;
22-
double beamStart = 500.; //[cm]
23-
double tokGauss = 1. / 0.1; // conversion from Tesla to kGauss
24-
25-
bool isMagAbs = true;
26-
27-
// ***********************
28-
// LAYOUT 1
29-
// ***********************
30-
3123
// RADIUS
32-
Rc = 170.; //[cm] — R_out_coil per Ian DetectorConstruction.cc; confirmed by A. Ortiz definition
33-
R1 = 220.; //[cm]
34-
R2 = 290.; //[cm]
24+
static constexpr double Rc = 170.; // [cm] — R_out_coil per Ian DetectorConstruction.cc; confirmed by A. Ortiz definition
25+
static constexpr double R1 = 220.; // [cm]
26+
static constexpr double R2 = 290.; // [cm]
3527

36-
// To set the B2
37-
B1 = 2.; //[T]
38-
B2 = -B1 * Rc * Rc / (R2 * R2 - R1 * R1); //[T] — B1 in numerator, confirmed by A. Ortiz Aug 2026
28+
// FIELD
29+
static constexpr double B1 = 2.; // [T]
30+
static constexpr double B2 = -B1 * Rc * Rc / (R2 * R2 - R1 * R1); // [T] — B1 in numerator, confirmed by A. Ortiz Aug 2026
31+
static constexpr double beamStart = 500.; // [cm]
32+
static constexpr double tokGauss = 1. / 0.1; // conversion from Tesla to kGauss
3933

40-
if ((abs(x[2]) <= beamStart) && (sqrt(x[0] * x[0] + x[1] * x[1]) < Rc)) {
34+
static constexpr bool isMagAbs = true;
35+
36+
const double r = sqrt(x[0] * x[0] + x[1] * x[1]);
37+
if ((abs(x[2]) <= beamStart) && (r < Rc)) { // We are inside of the central region
4138
b[0] = 0.;
4239
b[1] = 0.;
4340
b[2] = B1 * tokGauss;
44-
} else if ((abs(x[2]) <= beamStart) &&
45-
(sqrt(x[0] * x[0] + x[1] * x[1]) >= Rc &&
46-
sqrt(x[0] * x[0] + x[1] * x[1]) < R1)) {
41+
} else if ((abs(x[2]) <= beamStart) && (r >= Rc && r < R1)) { // We are in the transition region
4742
b[0] = 0.;
4843
b[1] = 0.;
4944
b[2] = 0.;
50-
} else if ((abs(x[2]) <= beamStart) &&
51-
(sqrt(x[0] * x[0] + x[1] * x[1]) >= R1 &&
52-
sqrt(x[0] * x[0] + x[1] * x[1]) < R2)) {
45+
} else if ((abs(x[2]) <= beamStart) && (r >= R1 && r < R2)) { // We are within the magnet
5346
b[0] = 0.;
5447
b[1] = 0.;
5548
if (isMagAbs) {
5649
b[2] = B2 * tokGauss;
5750
} else {
5851
b[2] = 0.;
5952
}
60-
} else {
53+
} else { // We are outside of the magnet
6154
b[0] = 0.;
6255
b[1] = 0.;
6356
b[2] = 0.;
6457
}
6558
};
6659
}
60+
61+
void ALICE3Field()
62+
{
63+
gStyle->SetPalette(kRainBow);
64+
gStyle->SetNumberContours(255);
65+
66+
auto fieldFunc = field();
67+
// RZ plane visualization
68+
TCanvas* cRZ = new TCanvas("cRZ", "Field in RZ plane", 800, 800);
69+
gPad->SetRightMargin(0.15);
70+
TH2F* hRZ = new TH2F("hRZ", "Magnetic Field B_z in RZ plane;Z [m];R [m];B_{z} [kGauss]", 100, -10, 10, 100, -5, 5);
71+
hRZ->SetBit(TH1::kNoStats); // disable stats box
72+
for (int i = 1; i <= hRZ->GetNbinsX(); i++) {
73+
const double Z = hRZ->GetXaxis()->GetBinCenter(i);
74+
for (int j = 1; j <= hRZ->GetNbinsY(); j++) {
75+
const double R = hRZ->GetYaxis()->GetBinCenter(j);
76+
const double pos[3] = {R * 100, 0, Z * 100}; // convert to cm
77+
double b[3] = {0, 0, 0};
78+
fieldFunc(pos, b);
79+
hRZ->SetBinContent(i, j, b[2]);
80+
}
81+
}
82+
83+
hRZ->GetZaxis()->SetRangeUser(-30, 30);
84+
hRZ->Draw("COLZ");
85+
cRZ->Update();
86+
87+
// XY plane visualization
88+
TCanvas* cXY = new TCanvas("cXY", "Field in XY plane", 800, 800);
89+
gPad->SetRightMargin(0.15);
90+
TH2F* hXY = new TH2F("hXY", "Magnetic Field B_z in XY plane;X [m];Y [m];B_{z} [kGauss]", 100, -5, 5, 100, -5, 5);
91+
hXY->SetBit(TH1::kNoStats); // disable stats box
92+
93+
for (int i = 1; i <= hXY->GetNbinsX(); i++) {
94+
const double X = hXY->GetXaxis()->GetBinCenter(i);
95+
for (int j = 1; j <= hXY->GetNbinsY(); j++) {
96+
const double Y = hXY->GetYaxis()->GetBinCenter(j);
97+
const double pos[3] = {X * 100, Y * 100, 0}; // convert to cm
98+
double b[3] = {0, 0, 0};
99+
fieldFunc(pos, b);
100+
hXY->SetBinContent(i, j, b[2]);
101+
}
102+
}
103+
104+
hXY->GetZaxis()->SetRangeUser(-30, 30);
105+
hXY->Draw("COLZ");
106+
cXY->Update();
107+
}

‎Detectors/Upgrades/ALICE3/macros/ALICE3FieldShortMagnet.C‎

Lines changed: 107 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,107 @@
1+
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
2+
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
3+
// All rights not expressly granted are reserved.
4+
//
5+
// This software is distributed under the terms of the GNU General Public
6+
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
7+
//
8+
// In applying this license CERN does not waive the privileges and immunities
9+
// granted to it by virtue of its status as an Intergovernmental Organization
10+
// or submit itself to any jurisdiction.
11+
//
12+
// Author: J. E. Munoz Mendez jesus.munoz@cern.ch
13+
14+
#include <functional>
15+
#include <cmath>
16+
#include <TCanvas.h>
17+
#include <TH2F.h>
18+
#include <TStyle.h>
19+
20+
std::function<void(const double*, double*)> field()
21+
{
22+
return [](const double* x, double* b) {
23+
// RADIUS
24+
static constexpr double Rc = 170.; // [cm] — R_out_coil per Ian DetectorConstruction.cc; confirmed by A. Ortiz definition
25+
static constexpr double R1 = 220.; // [cm]
26+
static constexpr double R2 = 290.; // [cm]
27+
28+
// FIELD
29+
static constexpr double B1 = 2.; // [T]
30+
static constexpr double B2 = -B1 * Rc * Rc / (R2 * R2 - R1 * R1); // [T] — B1 in numerator, confirmed by A. Ortiz Aug 2026
31+
static constexpr double beamStart = 370.; // [cm]
32+
static constexpr double tokGauss = 1. / 0.1; // conversion from Tesla to kGauss
33+
34+
static constexpr bool isMagAbs = true;
35+
36+
const double r = sqrt(x[0] * x[0] + x[1] * x[1]);
37+
if ((abs(x[2]) <= beamStart) && (r < Rc)) { // We are inside of the central region
38+
b[0] = 0.;
39+
b[1] = 0.;
40+
b[2] = B1 * tokGauss;
41+
} else if ((abs(x[2]) <= beamStart) && (r >= Rc && r < R1)) { // We are in the transition region
42+
b[0] = 0.;
43+
b[1] = 0.;
44+
b[2] = 0.;
45+
} else if ((abs(x[2]) <= beamStart) && (r >= R1 && r < R2)) { // We are within the magnet
46+
b[0] = 0.;
47+
b[1] = 0.;
48+
if (isMagAbs) {
49+
b[2] = B2 * tokGauss;
50+
} else {
51+
b[2] = 0.;
52+
}
53+
} else { // We are outside of the magnet
54+
b[0] = 0.;
55+
b[1] = 0.;
56+
b[2] = 0.;
57+
}
58+
};
59+
}
60+
61+
void ALICE3FieldShortMagnet()
62+
{
63+
gStyle->SetPalette(kRainBow);
64+
gStyle->SetNumberContours(255);
65+
66+
auto fieldFunc = field();
67+
// RZ plane visualization
68+
TCanvas* cRZ = new TCanvas("cRZ", "Field in RZ plane", 800, 800);
69+
gPad->SetRightMargin(0.15);
70+
TH2F* hRZ = new TH2F("hRZ", "Magnetic Field B_z in RZ plane;Z [m];R [m];B_{z} [kGauss]", 100, -10, 10, 100, -5, 5);
71+
hRZ->SetBit(TH1::kNoStats); // disable stats box
72+
for (int i = 1; i <= hRZ->GetNbinsX(); i++) {
73+
const double Z = hRZ->GetXaxis()->GetBinCenter(i);
74+
for (int j = 1; j <= hRZ->GetNbinsY(); j++) {
75+
const double R = hRZ->GetYaxis()->GetBinCenter(j);
76+
const double pos[3] = {R * 100, 0, Z * 100}; // convert to cm
77+
double b[3] = {0, 0, 0};
78+
fieldFunc(pos, b);
79+
hRZ->SetBinContent(i, j, b[2]);
80+
}
81+
}
82+
83+
hRZ->GetZaxis()->SetRangeUser(-30, 30);
84+
hRZ->Draw("COLZ");
85+
cRZ->Update();
86+
87+
// XY plane visualization
88+
TCanvas* cXY = new TCanvas("cXY", "Field in XY plane", 800, 800);
89+
gPad->SetRightMargin(0.15);
90+
TH2F* hXY = new TH2F("hXY", "Magnetic Field B_z in XY plane;X [m];Y [m];B_{z} [kGauss]", 100, -5, 5, 100, -5, 5);
91+
hXY->SetBit(TH1::kNoStats); // disable stats box
92+
93+
for (int i = 1; i <= hXY->GetNbinsX(); i++) {
94+
const double X = hXY->GetXaxis()->GetBinCenter(i);
95+
for (int j = 1; j <= hXY->GetNbinsY(); j++) {
96+
const double Y = hXY->GetYaxis()->GetBinCenter(j);
97+
const double pos[3] = {X * 100, Y * 100, 0}; // convert to cm
98+
double b[3] = {0, 0, 0};
99+
fieldFunc(pos, b);
100+
hXY->SetBinContent(i, j, b[2]);
101+
}
102+
}
103+
104+
hXY->GetZaxis()->SetRangeUser(-30, 30);
105+
hXY->Draw("COLZ");
106+
cXY->Update();
107+
}

‎Detectors/Upgrades/ALICE3/macros/CMakeLists.txt‎

Lines changed: 7 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -13,4 +13,10 @@ o2_add_test_root_macro(scanXX0.C
1313
LABELS alice3)
1414

1515
o2_add_test_root_macro(plotHits.C
16-
LABELS alice3)
16+
LABELS alice3)
17+
18+
o2_add_test_root_macro(ALICE3FieldShortMagnet.C
19+
LABELS alice3)
20+
21+
o2_add_test_root_macro(ALICE3Field.C
22+
LABELS alice3)

0 commit comments

Comments
 (0)