From defedd1aad5d2c3f39943da7138859ea0680117d Mon Sep 17 00:00:00 2001 From: olantwin Date: Wed, 20 Nov 2024 17:04:18 +0000 Subject: [PATCH] =?UTF-8?q?Deploying=20to=20gh-pages=20from=20=20@=2085dfb?= =?UTF-8?q?34206ef8edf46e2557c1817ba9e83234958=20=F0=9F=9A=80?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- ShipBFieldMap_8cxx_source.html | 769 +++++++++++++++++---------------- ShipBFieldMap_8h_source.html | 18 +- classShipBFieldMap.html | 769 +++++++++++++++++---------------- md_CHANGELOG.html | 1 + 4 files changed, 782 insertions(+), 775 deletions(-) diff --git a/ShipBFieldMap_8cxx_source.html b/ShipBFieldMap_8cxx_source.html index d0fc27c83..531047ed8 100644 --- a/ShipBFieldMap_8cxx_source.html +++ b/ShipBFieldMap_8cxx_source.html @@ -209,407 +209,410 @@
127  // 4. x < 0, y < 0: Bx = Bx
128 
129  Float_t BxSign(1.0);
-
130 
-
131  if (isSymmetric_) {
-
132 
-
133  // The field map co-ordinates only contain x > 0 and y > 0, i.e. we
-
134  // are using x-y quadrant symmetry. If the local x or y coordinates
-
135  // are negative we need to change their sign and keep track of the
-
136  // adjusted sign of Bx which we use as a multiplication factor at the end
-
137  if (x < 0.0) {
-
138  x = -x; BxSign *= -1.0;
-
139  }
-
140 
-
141  if (y < 0.0) {
-
142  y = -y; BxSign *= -1.0;
-
143  }
-
144 
-
145  }
-
146 
-
147  // Initialise the B field components to zero
-
148  B[0] = 0.0;
-
149  B[1] = 0.0;
-
150  B[2] = 0.0;
-
151 
-
152  // First check to see if we are inside the field map range
-
153  Bool_t inside = this->insideRange(x, y, z);
-
154  if (inside == kFALSE) {return;}
-
155 
-
156  // Find the neighbouring bins for the given point
-
157  binPair xBinInfo = this->getBinInfo(x, ShipBFieldMap::xAxis);
-
158  binPair yBinInfo = this->getBinInfo(y, ShipBFieldMap::yAxis);
-
159  binPair zBinInfo = this->getBinInfo(z, ShipBFieldMap::zAxis);
-
160 
-
161  Int_t iX = xBinInfo.first;
-
162  Int_t iY = yBinInfo.first;
-
163  Int_t iZ = zBinInfo.first;
-
164 
-
165  // Check that the bins are valid
-
166  if (iX == -1 || iY == -1 || iZ == -1) {return;}
+
130  Float_t BzSign = 1;
+
131 
+
132  if (isSymmetric_) {
+
133 
+
134  // The field map co-ordinates only contain x > 0 and y > 0, i.e. we
+
135  // are using x-y quadrant symmetry. If the local x or y coordinates
+
136  // are negative we need to change their sign and keep track of the
+
137  // adjusted sign of Bx which we use as a multiplication factor at the end
+
138  if (x < 0.0) {
+
139  x = -x; BxSign *= -1.0;
+
140  }
+
141 
+
142  if (y < 0.0) {
+
143  y = -y;
+
144  BxSign *= -1.0;
+
145  BzSign = -1.0;
+
146  }
+
147 
+
148  }
+
149 
+
150  // Initialise the B field components to zero
+
151  B[0] = 0.0;
+
152  B[1] = 0.0;
+
153  B[2] = 0.0;
+
154 
+
155  // First check to see if we are inside the field map range
+
156  Bool_t inside = this->insideRange(x, y, z);
+
157  if (inside == kFALSE) {return;}
+
158 
+
159  // Find the neighbouring bins for the given point
+
160  binPair xBinInfo = this->getBinInfo(x, ShipBFieldMap::xAxis);
+
161  binPair yBinInfo = this->getBinInfo(y, ShipBFieldMap::yAxis);
+
162  binPair zBinInfo = this->getBinInfo(z, ShipBFieldMap::zAxis);
+
163 
+
164  Int_t iX = xBinInfo.first;
+
165  Int_t iY = yBinInfo.first;
+
166  Int_t iZ = zBinInfo.first;
167 
-
168  // Get the various neighbouring bin entries
-
169  Int_t iX1(iX + 1);
-
170  Int_t iY1(iY + 1);
-
171  Int_t iZ1(iZ + 1);
-
172 
-
173  binA_ = this->getMapBin(iX, iY, iZ);
-
174  binB_ = this->getMapBin(iX1, iY, iZ);
-
175  binC_ = this->getMapBin(iX, iY1, iZ);
-
176  binD_ = this->getMapBin(iX1, iY1, iZ);
-
177  binE_ = this->getMapBin(iX, iY, iZ1);
-
178  binF_ = this->getMapBin(iX1, iY, iZ1);
-
179  binG_ = this->getMapBin(iX, iY1, iZ1);
-
180  binH_ = this->getMapBin(iX1, iY1, iZ1);
-
181 
-
182  // Retrieve the fractional bin distances
-
183  xFrac_ = xBinInfo.second;
-
184  yFrac_ = yBinInfo.second;
-
185  zFrac_ = zBinInfo.second;
-
186 
-
187  // Set the complimentary fractional bin distances
-
188  xFrac1_ = 1.0 - xFrac_;
-
189  yFrac1_ = 1.0 - yFrac_;
-
190  zFrac1_ = 1.0 - zFrac_;
-
191 
-
192  // Finally get the magnetic field components using trilinear interpolation
-
193  // and scale with the appropriate multiplication factor (default = 1.0)
-
194  B[0] = this->BInterCalc(ShipBFieldMap::xAxis)*scale_*BxSign;
-
195  B[1] = this->BInterCalc(ShipBFieldMap::yAxis)*scale_;
-
196  B[2] = this->BInterCalc(ShipBFieldMap::zAxis)*scale_;
-
197 
-
198 }
-
199 
-
200 void ShipBFieldMap::initialise()
-
201 {
+
168  // Check that the bins are valid
+
169  if (iX == -1 || iY == -1 || iZ == -1) {return;}
+
170 
+
171  // Get the various neighbouring bin entries
+
172  Int_t iX1(iX + 1);
+
173  Int_t iY1(iY + 1);
+
174  Int_t iZ1(iZ + 1);
+
175 
+
176  binA_ = this->getMapBin(iX, iY, iZ);
+
177  binB_ = this->getMapBin(iX1, iY, iZ);
+
178  binC_ = this->getMapBin(iX, iY1, iZ);
+
179  binD_ = this->getMapBin(iX1, iY1, iZ);
+
180  binE_ = this->getMapBin(iX, iY, iZ1);
+
181  binF_ = this->getMapBin(iX1, iY, iZ1);
+
182  binG_ = this->getMapBin(iX, iY1, iZ1);
+
183  binH_ = this->getMapBin(iX1, iY1, iZ1);
+
184 
+
185  // Retrieve the fractional bin distances
+
186  xFrac_ = xBinInfo.second;
+
187  yFrac_ = yBinInfo.second;
+
188  zFrac_ = zBinInfo.second;
+
189 
+
190  // Set the complimentary fractional bin distances
+
191  xFrac1_ = 1.0 - xFrac_;
+
192  yFrac1_ = 1.0 - yFrac_;
+
193  zFrac1_ = 1.0 - zFrac_;
+
194 
+
195  // Finally get the magnetic field components using trilinear interpolation
+
196  // and scale with the appropriate multiplication factor (default = 1.0)
+
197  B[0] = this->BInterCalc(ShipBFieldMap::xAxis)*scale_*BxSign;
+
198  B[1] = this->BInterCalc(ShipBFieldMap::yAxis)*scale_;
+
199  B[2] = this->BInterCalc(ShipBFieldMap::zAxis) * scale_ * BzSign;
+
200 
+
201 }
202 
-
203  if (initialised_ == kFALSE) {
-
204 
-
205  if (isCopy_ == kFALSE) {this->readMapFile();}
-
206 
-
207  // Set the global co-ordinate translation and rotation info
-
208  if (fabs(phi_) > 1e-6 || fabs(theta_) > 1e-6 || fabs(psi_) > 1e-6) {
+
203 void ShipBFieldMap::initialise()
+
204 {
+
205 
+
206  if (initialised_ == kFALSE) {
+
207 
+
208  if (isCopy_ == kFALSE) {this->readMapFile();}
209 
-
210  // We have non-zero rotation angles. Create a combined translation and rotation
-
211  TGeoTranslation tr("offsets", xOffset_, yOffset_, zOffset_);
-
212  TGeoRotation rot("angles", phi_, theta_, psi_);
-
213  theTrans_ = new TGeoCombiTrans(tr, rot);
-
214 
-
215  } else {
-
216 
-
217  // We only need a translation
-
218  theTrans_ = new TGeoTranslation("offsets", xOffset_, yOffset_, zOffset_);
+
210  // Set the global co-ordinate translation and rotation info
+
211  if (fabs(phi_) > 1e-6 || fabs(theta_) > 1e-6 || fabs(psi_) > 1e-6) {
+
212 
+
213  // We have non-zero rotation angles. Create a combined translation and rotation
+
214  TGeoTranslation tr("offsets", xOffset_, yOffset_, zOffset_);
+
215  TGeoRotation rot("angles", phi_, theta_, psi_);
+
216  theTrans_ = new TGeoCombiTrans(tr, rot);
+
217 
+
218  } else {
219 
-
220  }
-
221 
-
222  initialised_ = kTRUE;
-
223 
-
224  }
-
225 
-
226 }
-
227 
-
228 void ShipBFieldMap::readMapFile()
-
229 {
+
220  // We only need a translation
+
221  theTrans_ = new TGeoTranslation("offsets", xOffset_, yOffset_, zOffset_);
+
222 
+
223  }
+
224 
+
225  initialised_ = kTRUE;
+
226 
+
227  }
+
228 
+
229 }
230 
-
231  std::cout<<"ShipBFieldMap::readMapFile() creating field "<<this->GetName()
-
232  <<" using file "<<mapFileName_<<std::endl;
+
231 void ShipBFieldMap::readMapFile()
+
232 {
233 
-
234  // Check to see if we have a ROOT file
-
235  if (mapFileName_.find(".root") != std::string::npos) {
+
234  std::cout<<"ShipBFieldMap::readMapFile() creating field "<<this->GetName()
+
235  <<" using file "<<mapFileName_<<std::endl;
236 
-
237  this->readRootFile();
-
238 
-
239  } else {
-
240 
-
241  this->readTextFile();
-
242 
-
243  }
-
244 
-
245 }
-
246 
-
247 void ShipBFieldMap::readRootFile() {
-
248 
-
249  TFile* theFile = TFile::Open(mapFileName_.c_str());
-
250 
-
251  if (!theFile) {
-
252  std::cout<<"ShipBFieldMap: could not find the file "<<mapFileName_<<std::endl;
-
253  return;
-
254  }
-
255 
-
256  // Coordinate ranges
-
257  TTree* rTree = dynamic_cast<TTree*>(theFile->Get("Range"));
-
258  if (!rTree) {
-
259  std::cout<<"ShipBFieldMap: could not find Range tree in "<<mapFileName_<<std::endl;
-
260  return;
-
261  }
-
262 
-
263  rTree->SetBranchAddress("xMin", &xMin_);
-
264  rTree->SetBranchAddress("xMax", &xMax_);
-
265  rTree->SetBranchAddress("dx", &dx_);
-
266  rTree->SetBranchAddress("yMin", &yMin_);
-
267  rTree->SetBranchAddress("yMax", &yMax_);
-
268  rTree->SetBranchAddress("dy", &dy_);
-
269  rTree->SetBranchAddress("zMin", &zMin_);
-
270  rTree->SetBranchAddress("zMax", &zMax_);
-
271  rTree->SetBranchAddress("dz", &dz_);
-
272 
-
273  // Fill the ranges
-
274  rTree->GetEntry(0);
+
237  // Check to see if we have a ROOT file
+
238  if (mapFileName_.find(".root") != std::string::npos) {
+
239 
+
240  this->readRootFile();
+
241 
+
242  } else {
+
243 
+
244  this->readTextFile();
+
245 
+
246  }
+
247 
+
248 }
+
249 
+
250 void ShipBFieldMap::readRootFile() {
+
251 
+
252  TFile* theFile = TFile::Open(mapFileName_.c_str());
+
253 
+
254  if (!theFile) {
+
255  std::cout<<"ShipBFieldMap: could not find the file "<<mapFileName_<<std::endl;
+
256  return;
+
257  }
+
258 
+
259  // Coordinate ranges
+
260  TTree* rTree = dynamic_cast<TTree*>(theFile->Get("Range"));
+
261  if (!rTree) {
+
262  std::cout<<"ShipBFieldMap: could not find Range tree in "<<mapFileName_<<std::endl;
+
263  return;
+
264  }
+
265 
+
266  rTree->SetBranchAddress("xMin", &xMin_);
+
267  rTree->SetBranchAddress("xMax", &xMax_);
+
268  rTree->SetBranchAddress("dx", &dx_);
+
269  rTree->SetBranchAddress("yMin", &yMin_);
+
270  rTree->SetBranchAddress("yMax", &yMax_);
+
271  rTree->SetBranchAddress("dy", &dy_);
+
272  rTree->SetBranchAddress("zMin", &zMin_);
+
273  rTree->SetBranchAddress("zMax", &zMax_);
+
274  rTree->SetBranchAddress("dz", &dz_);
275 
-
276  this->setLimits();
-
277 
-
278  // Make sure we don't have a copy
-
279  if (isCopy_ == kFALSE) {
+
276  // Fill the ranges
+
277  rTree->GetEntry(0);
+
278 
+
279  this->setLimits();
280 
-
281  // The data is expected to contain Bx,By,Bz data values
-
282  // in ascending z,y,x co-ordinate order
+
281  // Make sure we don't have a copy
+
282  if (isCopy_ == kFALSE) {
283 
-
284  TTree* dTree = dynamic_cast<TTree*>(theFile->Get("Data"));
-
285  if (!dTree) {
-
286  std::cout<<"ShipBFieldMap: could not find Data tree in "<<mapFileName_<<std::endl;
-
287  return;
-
288  }
-
289 
-
290  Float_t Bx, By, Bz;
-
291  // Only enable the field components
-
292  dTree->SetBranchStatus("*", 0);
-
293  dTree->SetBranchStatus("Bx", 1);
-
294  dTree->SetBranchStatus("By", 1);
-
295  dTree->SetBranchStatus("Bz", 1);
-
296 
-
297  dTree->SetBranchAddress("Bx", &Bx);
-
298  dTree->SetBranchAddress("By", &By);
-
299  dTree->SetBranchAddress("Bz", &Bz);
-
300 
-
301  Int_t nEntries = dTree->GetEntries();
-
302  if (nEntries != N_) {
-
303  std::cout<<"Expected "<<N_<<" field map entries but found "<<nEntries<<std::endl;
-
304  nEntries = 0;
-
305  }
-
306 
-
307  fieldMap_->reserve(nEntries);
-
308 
-
309  for (Int_t i = 0; i < nEntries; i++) {
-
310 
-
311  dTree->GetEntry(i);
-
312 
-
313  // B field values are in Tesla. This means these values are multiplied
-
314  // by a factor of ten since both FairRoot and the VMC interface use kGauss
-
315  Bx *= Tesla_;
-
316  By *= Tesla_;
-
317  Bz *= Tesla_;
-
318 
-
319  // Store the B field 3-vector
-
320  std::vector<Float_t> BVector(3);
-
321  BVector[0] = Bx; BVector[1] = By; BVector[2] = Bz;
-
322  fieldMap_->push_back(BVector);
-
323 
-
324  }
-
325 
-
326  }
-
327 
-
328  theFile->Close();
-
329 
-
330 }
-
331 
-
332 void ShipBFieldMap::readTextFile() {
-
333 
-
334  std::ifstream getData(mapFileName_.c_str());
-
335 
-
336  std::string tmpString("");
-
337 
-
338  getData >> tmpString >> xMin_ >> xMax_ >> dx_
-
339  >> yMin_ >> yMax_ >> dy_ >> zMin_ >> zMax_ >> dz_;
+
284  // The data is expected to contain Bx,By,Bz data values
+
285  // in ascending z,y,x co-ordinate order
+
286 
+
287  TTree* dTree = dynamic_cast<TTree*>(theFile->Get("Data"));
+
288  if (!dTree) {
+
289  std::cout<<"ShipBFieldMap: could not find Data tree in "<<mapFileName_<<std::endl;
+
290  return;
+
291  }
+
292 
+
293  Float_t Bx, By, Bz;
+
294  // Only enable the field components
+
295  dTree->SetBranchStatus("*", 0);
+
296  dTree->SetBranchStatus("Bx", 1);
+
297  dTree->SetBranchStatus("By", 1);
+
298  dTree->SetBranchStatus("Bz", 1);
+
299 
+
300  dTree->SetBranchAddress("Bx", &Bx);
+
301  dTree->SetBranchAddress("By", &By);
+
302  dTree->SetBranchAddress("Bz", &Bz);
+
303 
+
304  Int_t nEntries = dTree->GetEntries();
+
305  if (nEntries != N_) {
+
306  std::cout<<"Expected "<<N_<<" field map entries but found "<<nEntries<<std::endl;
+
307  nEntries = 0;
+
308  }
+
309 
+
310  fieldMap_->reserve(nEntries);
+
311 
+
312  for (Int_t i = 0; i < nEntries; i++) {
+
313 
+
314  dTree->GetEntry(i);
+
315 
+
316  // B field values are in Tesla. This means these values are multiplied
+
317  // by a factor of ten since both FairRoot and the VMC interface use kGauss
+
318  Bx *= Tesla_;
+
319  By *= Tesla_;
+
320  Bz *= Tesla_;
+
321 
+
322  // Store the B field 3-vector
+
323  std::vector<Float_t> BVector(3);
+
324  BVector[0] = Bx; BVector[1] = By; BVector[2] = Bz;
+
325  fieldMap_->push_back(BVector);
+
326 
+
327  }
+
328 
+
329  }
+
330 
+
331  theFile->Close();
+
332 
+
333 }
+
334 
+
335 void ShipBFieldMap::readTextFile() {
+
336 
+
337  std::ifstream getData(mapFileName_.c_str());
+
338 
+
339  std::string tmpString("");
340 
-
341  this->setLimits();
-
342 
-
343  // Check to see if this object is a "copy"
-
344  if (isCopy_ == kFALSE) {
+
341  getData >> tmpString >> xMin_ >> xMax_ >> dx_
+
342  >> yMin_ >> yMax_ >> dy_ >> zMin_ >> zMax_ >> dz_;
+
343 
+
344  this->setLimits();
345 
-
346  // Read the field map and store the magnetic field components
-
347 
-
348  // Second line can be ignored since they are
-
349  // just labels for the data columns for readability
-
350  getData >> tmpString >> tmpString >> tmpString;
-
351 
-
352  // The remaining lines contain Bx,By,Bz data values
-
353  // in ascending z,y,x co-ord order
-
354  fieldMap_->reserve(N_);
-
355 
-
356  Float_t Bx(0.0), By(0.0), Bz(0.0);
-
357 
-
358  for (Int_t i = 0; i < N_; i++) {
-
359 
-
360  getData >> Bx >> By >> Bz;
-
361 
-
362  // B field values are in Tesla. This means these values are multiplied
-
363  // by a factor of ten since both FairRoot and the VMC interface use kGauss
-
364  Bx *= Tesla_;
-
365  By *= Tesla_;
-
366  Bz *= Tesla_;
-
367 
-
368  // Store the B field 3-vector
-
369  std::vector<Float_t> BVector(3);
-
370  BVector[0] = Bx; BVector[1] = By; BVector[2] = Bz;
-
371  fieldMap_->push_back(BVector);
-
372 
-
373  }
-
374 
-
375  }
-
376 
-
377  getData.close();
-
378 
-
379 }
-
380 
-
381 Bool_t ShipBFieldMap::insideRange(Float_t x, Float_t y, Float_t z)
-
382 {
+
346  // Check to see if this object is a "copy"
+
347  if (isCopy_ == kFALSE) {
+
348 
+
349  // Read the field map and store the magnetic field components
+
350 
+
351  // Second line can be ignored since they are
+
352  // just labels for the data columns for readability
+
353  getData >> tmpString >> tmpString >> tmpString;
+
354 
+
355  // The remaining lines contain Bx,By,Bz data values
+
356  // in ascending z,y,x co-ord order
+
357  fieldMap_->reserve(N_);
+
358 
+
359  Float_t Bx(0.0), By(0.0), Bz(0.0);
+
360 
+
361  for (Int_t i = 0; i < N_; i++) {
+
362 
+
363  getData >> Bx >> By >> Bz;
+
364 
+
365  // B field values are in Tesla. This means these values are multiplied
+
366  // by a factor of ten since both FairRoot and the VMC interface use kGauss
+
367  Bx *= Tesla_;
+
368  By *= Tesla_;
+
369  Bz *= Tesla_;
+
370 
+
371  // Store the B field 3-vector
+
372  std::vector<Float_t> BVector(3);
+
373  BVector[0] = Bx; BVector[1] = By; BVector[2] = Bz;
+
374  fieldMap_->push_back(BVector);
+
375 
+
376  }
+
377 
+
378  }
+
379 
+
380  getData.close();
+
381 
+
382 }
383 
-
384  Bool_t inside(kFALSE);
-
385 
-
386  if (x >= xMin_ && x <= xMax_ && y >= yMin_ &&
-
387  y <= yMax_ && z >= zMin_ && z <= zMax_) {inside = kTRUE;}
+
384 Bool_t ShipBFieldMap::insideRange(Float_t x, Float_t y, Float_t z)
+
385 {
+
386 
+
387  Bool_t inside(kFALSE);
388 
-
389  return inside;
-
390 
-
391 }
-
392 
-
393 void ShipBFieldMap::setLimits() {
-
394 
-
395  // Since the default SHIP distance unit is cm, we do not need to convert
-
396  // these map limits, i.e. cm = 1 already
+
389  if (x >= xMin_ && x <= xMax_ && y >= yMin_ &&
+
390  y <= yMax_ && z >= zMin_ && z <= zMax_) {inside = kTRUE;}
+
391 
+
392  return inside;
+
393 
+
394 }
+
395 
+
396 void ShipBFieldMap::setLimits() {
397 
-
398  xRange_ = xMax_ - xMin_;
-
399  yRange_ = yMax_ - yMin_;
-
400  zRange_ = zMax_ - zMin_;
-
401 
-
402  // Find the number of bins using the limits and bin sizes. The number of bins
-
403  // includes both the minimum and maximum values. To ensure correct rounding
-
404  // up to the nearest integer we need to add 1.5 not 1.0.
-
405  if (dx_ > 0.0) {
-
406  Nx_ = static_cast<Int_t>(((xMax_ - xMin_)/dx_) + 1.5);
-
407  }
-
408  if (dy_ > 0.0) {
-
409  Ny_ = static_cast<Int_t>(((yMax_ - yMin_)/dy_) + 1.5);
+
398  // Since the default SHIP distance unit is cm, we do not need to convert
+
399  // these map limits, i.e. cm = 1 already
+
400 
+
401  xRange_ = xMax_ - xMin_;
+
402  yRange_ = yMax_ - yMin_;
+
403  zRange_ = zMax_ - zMin_;
+
404 
+
405  // Find the number of bins using the limits and bin sizes. The number of bins
+
406  // includes both the minimum and maximum values. To ensure correct rounding
+
407  // up to the nearest integer we need to add 1.5 not 1.0.
+
408  if (dx_ > 0.0) {
+
409  Nx_ = static_cast<Int_t>(((xMax_ - xMin_)/dx_) + 1.5);
410  }
-
411  if (dz_ > 0.0) {
-
412  Nz_ = static_cast<Int_t>(((zMax_ - zMin_)/dz_) + 1.5);
+
411  if (dy_ > 0.0) {
+
412  Ny_ = static_cast<Int_t>(((yMax_ - yMin_)/dy_) + 1.5);
413  }
-
414 
-
415  N_ = Nx_*Ny_*Nz_;
-
416 
-
417  std::cout<<"x limits: "<<xMin_<<", "<<xMax_<<", dx = "<<dx_<<std::endl;
-
418  std::cout<<"y limits: "<<yMin_<<", "<<yMax_<<", dy = "<<dy_<<std::endl;
-
419  std::cout<<"z limits: "<<zMin_<<", "<<zMax_<<", dz = "<<dz_<<std::endl;
-
420 
-
421  std::cout<<"Offsets: x = "<<xOffset_<<", y = "<<yOffset_<<", z = "<<zOffset_<<std::endl;
-
422  std::cout<<"Angles : phi = "<<phi_<<", theta = "<<theta_<<", psi = "<<psi_<<std::endl;
+
414  if (dz_ > 0.0) {
+
415  Nz_ = static_cast<Int_t>(((zMax_ - zMin_)/dz_) + 1.5);
+
416  }
+
417 
+
418  N_ = Nx_*Ny_*Nz_;
+
419 
+
420  std::cout<<"x limits: "<<xMin_<<", "<<xMax_<<", dx = "<<dx_<<std::endl;
+
421  std::cout<<"y limits: "<<yMin_<<", "<<yMax_<<", dy = "<<dy_<<std::endl;
+
422  std::cout<<"z limits: "<<zMin_<<", "<<zMax_<<", dz = "<<dz_<<std::endl;
423 
-
424  std::cout<<"Total number of bins = "<<N_
-
425  <<"; Nx = "<<Nx_<<", Ny = "<<Ny_<<", Nz = "<<Nz_<<std::endl;
+
424  std::cout<<"Offsets: x = "<<xOffset_<<", y = "<<yOffset_<<", z = "<<zOffset_<<std::endl;
+
425  std::cout<<"Angles : phi = "<<phi_<<", theta = "<<theta_<<", psi = "<<psi_<<std::endl;
426 
-
427 }
-
428 
+
427  std::cout<<"Total number of bins = "<<N_
+
428  <<"; Nx = "<<Nx_<<", Ny = "<<Ny_<<", Nz = "<<Nz_<<std::endl;
429 
-
430 ShipBFieldMap::binPair ShipBFieldMap::getBinInfo(Float_t u, ShipBFieldMap::CoordAxis theAxis)
-
431 {
+
430 }
+
431 
432 
-
433  Float_t du(0.0), uMin(0.0), Nu(0);
-
434 
-
435  if (theAxis == ShipBFieldMap::xAxis) {
-
436  du = dx_; uMin = xMin_; Nu = Nx_;
-
437  } else if (theAxis == ShipBFieldMap::yAxis) {
-
438  du = dy_; uMin = yMin_; Nu = Ny_;
-
439  } else if (theAxis == ShipBFieldMap::zAxis) {
-
440  du = dz_; uMin = zMin_; Nu = Nz_;
-
441  }
-
442 
-
443  Int_t iBin(-1);
-
444  Float_t fracL(0.0);
+
433 ShipBFieldMap::binPair ShipBFieldMap::getBinInfo(Float_t u, ShipBFieldMap::CoordAxis theAxis)
+
434 {
+
435 
+
436  Float_t du(0.0), uMin(0.0), Nu(0);
+
437 
+
438  if (theAxis == ShipBFieldMap::xAxis) {
+
439  du = dx_; uMin = xMin_; Nu = Nx_;
+
440  } else if (theAxis == ShipBFieldMap::yAxis) {
+
441  du = dy_; uMin = yMin_; Nu = Ny_;
+
442  } else if (theAxis == ShipBFieldMap::zAxis) {
+
443  du = dz_; uMin = zMin_; Nu = Nz_;
+
444  }
445 
-
446  if (du > 1e-10) {
-
447 
-
448  // Get the number of fractional bin widths the point is from the first volume bin
-
449  Float_t dist = (u - uMin)/du;
-
450  // Get the integer equivalent of this distance, which is the bin number
-
451  iBin = static_cast<Int_t>(dist);
-
452  // Get the actual fractional distance of the point from the leftmost bin edge
-
453  fracL = (dist - iBin*1.0);
-
454 
-
455  }
-
456 
-
457  // Check that the bin number is valid
-
458  if (iBin < 0 || iBin >= Nu) {
-
459  iBin = -1; fracL = 0.0;
-
460  }
-
461 
-
462  // Return the bin number and fractional distance of the point from the leftmost bin edge
-
463  binPair binInfo(iBin, fracL);
+
446  Int_t iBin(-1);
+
447  Float_t fracL(0.0);
+
448 
+
449  if (du > 1e-10) {
+
450 
+
451  // Get the number of fractional bin widths the point is from the first volume bin
+
452  Float_t dist = (u - uMin)/du;
+
453  // Get the integer equivalent of this distance, which is the bin number
+
454  iBin = static_cast<Int_t>(dist);
+
455  // Get the actual fractional distance of the point from the leftmost bin edge
+
456  fracL = (dist - iBin*1.0);
+
457 
+
458  }
+
459 
+
460  // Check that the bin number is valid
+
461  if (iBin < 0 || iBin >= Nu) {
+
462  iBin = -1; fracL = 0.0;
+
463  }
464 
-
465  return binInfo;
-
466 
-
467 }
-
468 
-
469 Int_t ShipBFieldMap::getMapBin(Int_t iX, Int_t iY, Int_t iZ)
-
470 {
+
465  // Return the bin number and fractional distance of the point from the leftmost bin edge
+
466  binPair binInfo(iBin, fracL);
+
467 
+
468  return binInfo;
+
469 
+
470 }
471 
-
472  // Get the index of the map entry corresponding to the x,y,z bins.
-
473  // Remember that the map is ordered in ascending z, y, then x
+
472 Int_t ShipBFieldMap::getMapBin(Int_t iX, Int_t iY, Int_t iZ)
+
473 {
474 
-
475  Int_t index = (iX*Ny_ + iY)*Nz_ + iZ;
-
476  if (index < 0) {
-
477  index = 0;
-
478  } else if (index >= N_) {
-
479  index = N_-1;
-
480  }
-
481 
-
482  return index;
-
483 
-
484 }
-
485 
-
486 Float_t ShipBFieldMap::BInterCalc(CoordAxis theAxis)
-
487 {
+
475  // Get the index of the map entry corresponding to the x,y,z bins.
+
476  // Remember that the map is ordered in ascending z, y, then x
+
477 
+
478  Int_t index = (iX*Ny_ + iY)*Nz_ + iZ;
+
479  if (index < 0) {
+
480  index = 0;
+
481  } else if (index >= N_) {
+
482  index = N_-1;
+
483  }
+
484 
+
485  return index;
+
486 
+
487 }
488 
-
489  // Find the magnetic field component along theAxis using trilinear
-
490  // interpolation based on the current position and neighbouring bins
-
491  Float_t result(0.0);
-
492 
-
493  Int_t iAxis(0);
-
494 
-
495  if (theAxis == CoordAxis::yAxis) {
-
496  iAxis = 1;
-
497  } else if (theAxis == CoordAxis::zAxis) {
-
498  iAxis = 2;
-
499  }
-
500 
-
501  if (fieldMap_) {
-
502 
-
503  // Get the field component values for the neighbouring bins
-
504  Float_t A = (*fieldMap_)[binA_][iAxis];
-
505  Float_t B = (*fieldMap_)[binB_][iAxis];
-
506  Float_t C = (*fieldMap_)[binC_][iAxis];
-
507  Float_t D = (*fieldMap_)[binD_][iAxis];
-
508  Float_t E = (*fieldMap_)[binE_][iAxis];
-
509  Float_t F = (*fieldMap_)[binF_][iAxis];
-
510  Float_t G = (*fieldMap_)[binG_][iAxis];
-
511  Float_t H = (*fieldMap_)[binH_][iAxis];
-
512 
-
513  // Perform linear interpolation along x
-
514  Float_t F00 = A*xFrac1_ + B*xFrac_;
-
515  Float_t F10 = C*xFrac1_ + D*xFrac_;
-
516  Float_t F01 = E*xFrac1_ + F*xFrac_;
-
517  Float_t F11 = G*xFrac1_ + H*xFrac_;
-
518 
-
519  // Linear interpolation along y
-
520  Float_t F0 = F00*yFrac1_ + F10*yFrac_;
-
521  Float_t F1 = F01*yFrac1_ + F11*yFrac_;
-
522 
-
523  // Linear interpolation along z
-
524  result = F0*zFrac1_ + F1*zFrac_;
+
489 Float_t ShipBFieldMap::BInterCalc(CoordAxis theAxis)
+
490 {
+
491 
+
492  // Find the magnetic field component along theAxis using trilinear
+
493  // interpolation based on the current position and neighbouring bins
+
494  Float_t result(0.0);
+
495 
+
496  Int_t iAxis(0);
+
497 
+
498  if (theAxis == CoordAxis::yAxis) {
+
499  iAxis = 1;
+
500  } else if (theAxis == CoordAxis::zAxis) {
+
501  iAxis = 2;
+
502  }
+
503 
+
504  if (fieldMap_) {
+
505 
+
506  // Get the field component values for the neighbouring bins
+
507  Float_t A = (*fieldMap_)[binA_][iAxis];
+
508  Float_t B = (*fieldMap_)[binB_][iAxis];
+
509  Float_t C = (*fieldMap_)[binC_][iAxis];
+
510  Float_t D = (*fieldMap_)[binD_][iAxis];
+
511  Float_t E = (*fieldMap_)[binE_][iAxis];
+
512  Float_t F = (*fieldMap_)[binF_][iAxis];
+
513  Float_t G = (*fieldMap_)[binG_][iAxis];
+
514  Float_t H = (*fieldMap_)[binH_][iAxis];
+
515 
+
516  // Perform linear interpolation along x
+
517  Float_t F00 = A*xFrac1_ + B*xFrac_;
+
518  Float_t F10 = C*xFrac1_ + D*xFrac_;
+
519  Float_t F01 = E*xFrac1_ + F*xFrac_;
+
520  Float_t F11 = G*xFrac1_ + H*xFrac_;
+
521 
+
522  // Linear interpolation along y
+
523  Float_t F0 = F00*yFrac1_ + F10*yFrac_;
+
524  Float_t F1 = F01*yFrac1_ + F11*yFrac_;
525 
-
526  }
-
527 
-
528  return result;
-
529 
-
530 }
+
526  // Linear interpolation along z
+
527  result = F0*zFrac1_ + F1*zFrac_;
+
528 
+
529  }
+
530 
+
531  return result;
+
532 
+
533 }
Class that defines a (3d) magnetic field map (distances in cm, fields in tesla)
Definition: ShipBFieldMap.h:17
Float_t zFrac_
Fractional bin distance along z.
@@ -624,7 +627,7 @@
Float_t zFrac1_
Complimentary fractional bin distance along z.
ShipBFieldMap(const std::string &label, const std::string &mapFileName, Float_t xOffset=0.0, Float_t yOffset=0.0, Float_t zOffset=0.0, Float_t phi=0.0, Float_t theta=0.0, Float_t psi=0.0, Float_t scale=1.0, Bool_t isSymmetric=kFALSE)
Constructor.
Int_t binE_
Bin E for the trilinear interpolation.
-
Float_t BInterCalc(CoordAxis theAxis)
+
Float_t BInterCalc(CoordAxis theAxis)
Int_t Nz_
The number of bins along z.
Float_t Tesla_
Double converting Tesla to kiloGauss (for VMC/FairRoot B field units)
Int_t N_
The total number of bins.
@@ -634,7 +637,7 @@
Float_t yMin_
The minimum value of y for the map.
Float_t zOffset_
The z value of the positional offset for the map.
std::string mapFileName_
The name of the map file.
-
void setLimits()
+
void setLimits()
Bool_t initialised_
Initialisation boolean.
Int_t Ny_
The number of bins along y.
Float_t theta_
The second Euler rotation angle about the new x axis.
@@ -645,7 +648,7 @@
Int_t binG_
Bin G for the trilinear interpolation.
Int_t binH_
Bin H for the trilinear interpolation.
Float_t xMin_
The minimum value of x for the map.
-
void initialise()
Initialisation.
+
void initialise()
Initialisation.
Int_t binB_
Bin B for the trilinear interpolation.
Int_t binF_
Bin F for the trilinear interpolation.
CoordAxis
Enumeration to specify the co-ordinate type.
@@ -655,17 +658,17 @@
Int_t Nx_
The number of bins along x.
Float_t xRange_
The co-ordinate range along x.
Float_t yOffset_
The y value of the positional offset for the map.
-
void readTextFile()
Process the text file containing the field map data.
+
void readTextFile()
Process the text file containing the field map data.
Float_t scale_
The B field magnitude scaling factor.
TGeoMatrix * theTrans_
The combined translation and rotation transformation.
-
binPair getBinInfo(Float_t x, CoordAxis theAxis)
Get the bin number and fractional distance from the leftmost bin edge.
-
Bool_t insideRange(Float_t x, Float_t y, Float_t z)
Check to see if a point is within the map validity range.
+
binPair getBinInfo(Float_t x, CoordAxis theAxis)
Get the bin number and fractional distance from the leftmost bin edge.
+
Bool_t insideRange(Float_t x, Float_t y, Float_t z)
Check to see if a point is within the map validity range.
virtual ~ShipBFieldMap()
Destructor.
-
void readRootFile()
Process the ROOT file containing the field map data.
+
void readRootFile()
Process the ROOT file containing the field map data.
Float_t xFrac1_
Complimentary fractional bin distance along x.
Int_t binA_
Bin A for the trilinear interpolation.
-
void readMapFile()
Read the field map data and store the information internally.
-
Int_t getMapBin(Int_t iX, Int_t iY, Int_t iZ)
Find the vector entry of the field map data given the bins iX, iY and iZ.
+
void readMapFile()
Read the field map data and store the information internally.
+
Int_t getMapBin(Int_t iX, Int_t iY, Int_t iZ)
Find the vector entry of the field map data given the bins iX, iY and iZ.
Float_t yMax_
The maximum value of y for the map.
Float_t yRange_
The co-ordinate range along y.
Float_t zMax_
The maximum value of z for the map.
diff --git a/ShipBFieldMap_8h_source.html b/ShipBFieldMap_8h_source.html index e088407cd..900569a29 100644 --- a/ShipBFieldMap_8h_source.html +++ b/ShipBFieldMap_8h_source.html @@ -375,7 +375,7 @@
Float_t zFrac1_
Complimentary fractional bin distance along z.
ShipBFieldMap(const std::string &label, const std::string &mapFileName, Float_t xOffset=0.0, Float_t yOffset=0.0, Float_t zOffset=0.0, Float_t phi=0.0, Float_t theta=0.0, Float_t psi=0.0, Float_t scale=1.0, Bool_t isSymmetric=kFALSE)
Constructor.
Int_t binE_
Bin E for the trilinear interpolation.
-
Float_t BInterCalc(CoordAxis theAxis)
+
Float_t BInterCalc(CoordAxis theAxis)
Int_t Nz_
The number of bins along z.
Float_t Tesla_
Double converting Tesla to kiloGauss (for VMC/FairRoot B field units)
Float_t GetXRange() const
Get the x co-ordinate range for the map.
@@ -388,7 +388,7 @@
Float_t yMin_
The minimum value of y for the map.
Float_t zOffset_
The z value of the positional offset for the map.
std::string mapFileName_
The name of the map file.
-
void setLimits()
+
void setLimits()
Bool_t initialised_
Initialisation boolean.
void SetPhi(Float_t phi)
Set the first Euler rotation angle phi about the z axis.
Definition: ShipBFieldMap.h:99
Int_t Ny_
The number of bins along y.
@@ -403,7 +403,7 @@
Int_t binG_
Bin G for the trilinear interpolation.
Int_t binH_
Bin H for the trilinear interpolation.
Float_t xMin_
The minimum value of x for the map.
-
void initialise()
Initialisation.
+
void initialise()
Initialisation.
void SetYOffset(Float_t yValue)
Set the y global co-ordinate shift.
Definition: ShipBFieldMap.h:87
Int_t GetNz() const
Get the number of bins along z.
Int_t binB_
Bin B for the trilinear interpolation.
@@ -425,21 +425,21 @@
Float_t xRange_
The co-ordinate range along x.
Float_t yOffset_
The y value of the positional offset for the map.
ShipBFieldMap & operator=(const ShipBFieldMap &rhs)
Copy assignment operator not implemented.
-
void readTextFile()
Process the text file containing the field map data.
+
void readTextFile()
Process the text file containing the field map data.
Float_t scale_
The B field magnitude scaling factor.
TGeoMatrix * theTrans_
The combined translation and rotation transformation.
-
binPair getBinInfo(Float_t x, CoordAxis theAxis)
Get the bin number and fractional distance from the leftmost bin edge.
+
binPair getBinInfo(Float_t x, CoordAxis theAxis)
Get the bin number and fractional distance from the leftmost bin edge.
Float_t GetPhi() const
Get the first Euler rotation angle about the z axis for global positioning.
-
Bool_t insideRange(Float_t x, Float_t y, Float_t z)
Check to see if a point is within the map validity range.
+
Bool_t insideRange(Float_t x, Float_t y, Float_t z)
Check to see if a point is within the map validity range.
void SetTheta(Float_t theta)
Set the second Euler rotation angle theta about the new x axis.
virtual ~ShipBFieldMap()
Destructor.
-
void readRootFile()
Process the ROOT file containing the field map data.
+
void readRootFile()
Process the ROOT file containing the field map data.
Float_t GetdZ() const
Get the bin width along z for the map.
Float_t xFrac1_
Complimentary fractional bin distance along x.
Int_t binA_
Bin A for the trilinear interpolation.
-
void readMapFile()
Read the field map data and store the information internally.
+
void readMapFile()
Read the field map data and store the information internally.
Float_t GetXOffset() const
Get the x offset co-ordinate of the map for global positioning.
-
Int_t getMapBin(Int_t iX, Int_t iY, Int_t iZ)
Find the vector entry of the field map data given the bins iX, iY and iZ.
+
Int_t getMapBin(Int_t iX, Int_t iY, Int_t iZ)
Find the vector entry of the field map data given the bins iX, iY and iZ.
Float_t yMax_
The maximum value of y for the map.
Float_t yRange_
The co-ordinate range along y.
Float_t zMax_
The maximum value of z for the map.
diff --git a/classShipBFieldMap.html b/classShipBFieldMap.html index e97fd0e73..b0fca281b 100644 --- a/classShipBFieldMap.html +++ b/classShipBFieldMap.html @@ -649,7 +649,7 @@

Float_t theta_
The second Euler rotation angle about the new x axis.
Bool_t isCopy_
Flag to specify if we are a copy of the field map (just a change of global offsets)
Float_t xMin_
The minimum value of x for the map.
-
void initialise()
Initialisation.
+
void initialise()
Initialisation.
Int_t Nx_
The number of bins along x.
Float_t xRange_
The co-ordinate range along x.
Float_t yOffset_
The y value of the positional offset for the map.
@@ -895,51 +895,51 @@

Returns
the magnetic field component for the given axis
-

Definition at line 486 of file ShipBFieldMap.cxx.

-
487 {
-
488 
-
489  // Find the magnetic field component along theAxis using trilinear
-
490  // interpolation based on the current position and neighbouring bins
-
491  Float_t result(0.0);
-
492 
-
493  Int_t iAxis(0);
-
494 
-
495  if (theAxis == CoordAxis::yAxis) {
-
496  iAxis = 1;
-
497  } else if (theAxis == CoordAxis::zAxis) {
-
498  iAxis = 2;
-
499  }
-
500 
-
501  if (fieldMap_) {
-
502 
-
503  // Get the field component values for the neighbouring bins
-
504  Float_t A = (*fieldMap_)[binA_][iAxis];
-
505  Float_t B = (*fieldMap_)[binB_][iAxis];
-
506  Float_t C = (*fieldMap_)[binC_][iAxis];
-
507  Float_t D = (*fieldMap_)[binD_][iAxis];
-
508  Float_t E = (*fieldMap_)[binE_][iAxis];
-
509  Float_t F = (*fieldMap_)[binF_][iAxis];
-
510  Float_t G = (*fieldMap_)[binG_][iAxis];
-
511  Float_t H = (*fieldMap_)[binH_][iAxis];
-
512 
-
513  // Perform linear interpolation along x
-
514  Float_t F00 = A*xFrac1_ + B*xFrac_;
-
515  Float_t F10 = C*xFrac1_ + D*xFrac_;
-
516  Float_t F01 = E*xFrac1_ + F*xFrac_;
-
517  Float_t F11 = G*xFrac1_ + H*xFrac_;
-
518 
-
519  // Linear interpolation along y
-
520  Float_t F0 = F00*yFrac1_ + F10*yFrac_;
-
521  Float_t F1 = F01*yFrac1_ + F11*yFrac_;
-
522 
-
523  // Linear interpolation along z
-
524  result = F0*zFrac1_ + F1*zFrac_;
+

Definition at line 489 of file ShipBFieldMap.cxx.

+
490 {
+
491 
+
492  // Find the magnetic field component along theAxis using trilinear
+
493  // interpolation based on the current position and neighbouring bins
+
494  Float_t result(0.0);
+
495 
+
496  Int_t iAxis(0);
+
497 
+
498  if (theAxis == CoordAxis::yAxis) {
+
499  iAxis = 1;
+
500  } else if (theAxis == CoordAxis::zAxis) {
+
501  iAxis = 2;
+
502  }
+
503 
+
504  if (fieldMap_) {
+
505 
+
506  // Get the field component values for the neighbouring bins
+
507  Float_t A = (*fieldMap_)[binA_][iAxis];
+
508  Float_t B = (*fieldMap_)[binB_][iAxis];
+
509  Float_t C = (*fieldMap_)[binC_][iAxis];
+
510  Float_t D = (*fieldMap_)[binD_][iAxis];
+
511  Float_t E = (*fieldMap_)[binE_][iAxis];
+
512  Float_t F = (*fieldMap_)[binF_][iAxis];
+
513  Float_t G = (*fieldMap_)[binG_][iAxis];
+
514  Float_t H = (*fieldMap_)[binH_][iAxis];
+
515 
+
516  // Perform linear interpolation along x
+
517  Float_t F00 = A*xFrac1_ + B*xFrac_;
+
518  Float_t F10 = C*xFrac1_ + D*xFrac_;
+
519  Float_t F01 = E*xFrac1_ + F*xFrac_;
+
520  Float_t F11 = G*xFrac1_ + H*xFrac_;
+
521 
+
522  // Linear interpolation along y
+
523  Float_t F0 = F00*yFrac1_ + F10*yFrac_;
+
524  Float_t F1 = F01*yFrac1_ + F11*yFrac_;
525 
-
526  }
-
527 
-
528  return result;
-
529 
-
530 }
+
526  // Linear interpolation along z
+
527  result = F0*zFrac1_ + F1*zFrac_;
+
528 
+
529  }
+
530 
+
531  return result;
+
532 
+
533 }
Float_t zFrac_
Fractional bin distance along z.
Int_t binD_
Bin D for the trilinear interpolation.
Float_t zFrac1_
Complimentary fractional bin distance along z.
@@ -1062,80 +1062,83 @@

127  // 4. x < 0, y < 0: Bx = Bx

128 
129  Float_t BxSign(1.0);
-
130 
-
131  if (isSymmetric_) {
-
132 
-
133  // The field map co-ordinates only contain x > 0 and y > 0, i.e. we
-
134  // are using x-y quadrant symmetry. If the local x or y coordinates
-
135  // are negative we need to change their sign and keep track of the
-
136  // adjusted sign of Bx which we use as a multiplication factor at the end
-
137  if (x < 0.0) {
-
138  x = -x; BxSign *= -1.0;
-
139  }
-
140 
-
141  if (y < 0.0) {
-
142  y = -y; BxSign *= -1.0;
-
143  }
-
144 
-
145  }
-
146 
-
147  // Initialise the B field components to zero
-
148  B[0] = 0.0;
-
149  B[1] = 0.0;
-
150  B[2] = 0.0;
-
151 
-
152  // First check to see if we are inside the field map range
-
153  Bool_t inside = this->insideRange(x, y, z);
-
154  if (inside == kFALSE) {return;}
-
155 
-
156  // Find the neighbouring bins for the given point
-
157  binPair xBinInfo = this->getBinInfo(x, ShipBFieldMap::xAxis);
-
158  binPair yBinInfo = this->getBinInfo(y, ShipBFieldMap::yAxis);
-
159  binPair zBinInfo = this->getBinInfo(z, ShipBFieldMap::zAxis);
-
160 
-
161  Int_t iX = xBinInfo.first;
-
162  Int_t iY = yBinInfo.first;
-
163  Int_t iZ = zBinInfo.first;
-
164 
-
165  // Check that the bins are valid
-
166  if (iX == -1 || iY == -1 || iZ == -1) {return;}
+
130  Float_t BzSign = 1;
+
131 
+
132  if (isSymmetric_) {
+
133 
+
134  // The field map co-ordinates only contain x > 0 and y > 0, i.e. we
+
135  // are using x-y quadrant symmetry. If the local x or y coordinates
+
136  // are negative we need to change their sign and keep track of the
+
137  // adjusted sign of Bx which we use as a multiplication factor at the end
+
138  if (x < 0.0) {
+
139  x = -x; BxSign *= -1.0;
+
140  }
+
141 
+
142  if (y < 0.0) {
+
143  y = -y;
+
144  BxSign *= -1.0;
+
145  BzSign = -1.0;
+
146  }
+
147 
+
148  }
+
149 
+
150  // Initialise the B field components to zero
+
151  B[0] = 0.0;
+
152  B[1] = 0.0;
+
153  B[2] = 0.0;
+
154 
+
155  // First check to see if we are inside the field map range
+
156  Bool_t inside = this->insideRange(x, y, z);
+
157  if (inside == kFALSE) {return;}
+
158 
+
159  // Find the neighbouring bins for the given point
+
160  binPair xBinInfo = this->getBinInfo(x, ShipBFieldMap::xAxis);
+
161  binPair yBinInfo = this->getBinInfo(y, ShipBFieldMap::yAxis);
+
162  binPair zBinInfo = this->getBinInfo(z, ShipBFieldMap::zAxis);
+
163 
+
164  Int_t iX = xBinInfo.first;
+
165  Int_t iY = yBinInfo.first;
+
166  Int_t iZ = zBinInfo.first;
167 
-
168  // Get the various neighbouring bin entries
-
169  Int_t iX1(iX + 1);
-
170  Int_t iY1(iY + 1);
-
171  Int_t iZ1(iZ + 1);
-
172 
-
173  binA_ = this->getMapBin(iX, iY, iZ);
-
174  binB_ = this->getMapBin(iX1, iY, iZ);
-
175  binC_ = this->getMapBin(iX, iY1, iZ);
-
176  binD_ = this->getMapBin(iX1, iY1, iZ);
-
177  binE_ = this->getMapBin(iX, iY, iZ1);
-
178  binF_ = this->getMapBin(iX1, iY, iZ1);
-
179  binG_ = this->getMapBin(iX, iY1, iZ1);
-
180  binH_ = this->getMapBin(iX1, iY1, iZ1);
-
181 
-
182  // Retrieve the fractional bin distances
-
183  xFrac_ = xBinInfo.second;
-
184  yFrac_ = yBinInfo.second;
-
185  zFrac_ = zBinInfo.second;
-
186 
-
187  // Set the complimentary fractional bin distances
-
188  xFrac1_ = 1.0 - xFrac_;
-
189  yFrac1_ = 1.0 - yFrac_;
-
190  zFrac1_ = 1.0 - zFrac_;
-
191 
-
192  // Finally get the magnetic field components using trilinear interpolation
-
193  // and scale with the appropriate multiplication factor (default = 1.0)
-
194  B[0] = this->BInterCalc(ShipBFieldMap::xAxis)*scale_*BxSign;
- - -
197 
-
198 }
+
168  // Check that the bins are valid
+
169  if (iX == -1 || iY == -1 || iZ == -1) {return;}
+
170 
+
171  // Get the various neighbouring bin entries
+
172  Int_t iX1(iX + 1);
+
173  Int_t iY1(iY + 1);
+
174  Int_t iZ1(iZ + 1);
+
175 
+
176  binA_ = this->getMapBin(iX, iY, iZ);
+
177  binB_ = this->getMapBin(iX1, iY, iZ);
+
178  binC_ = this->getMapBin(iX, iY1, iZ);
+
179  binD_ = this->getMapBin(iX1, iY1, iZ);
+
180  binE_ = this->getMapBin(iX, iY, iZ1);
+
181  binF_ = this->getMapBin(iX1, iY, iZ1);
+
182  binG_ = this->getMapBin(iX, iY1, iZ1);
+
183  binH_ = this->getMapBin(iX1, iY1, iZ1);
+
184 
+
185  // Retrieve the fractional bin distances
+
186  xFrac_ = xBinInfo.second;
+
187  yFrac_ = yBinInfo.second;
+
188  zFrac_ = zBinInfo.second;
+
189 
+
190  // Set the complimentary fractional bin distances
+
191  xFrac1_ = 1.0 - xFrac_;
+
192  yFrac1_ = 1.0 - yFrac_;
+
193  zFrac1_ = 1.0 - zFrac_;
+
194 
+
195  // Finally get the magnetic field components using trilinear interpolation
+
196  // and scale with the appropriate multiplication factor (default = 1.0)
+
197  B[0] = this->BInterCalc(ShipBFieldMap::xAxis)*scale_*BxSign;
+ +
199  B[2] = this->BInterCalc(ShipBFieldMap::zAxis) * scale_ * BzSign;
+
200 
+
201 }
std::pair< Int_t, Float_t > binPair
Typedef for an int-double pair.
-
Float_t BInterCalc(CoordAxis theAxis)
-
binPair getBinInfo(Float_t x, CoordAxis theAxis)
Get the bin number and fractional distance from the leftmost bin edge.
-
Bool_t insideRange(Float_t x, Float_t y, Float_t z)
Check to see if a point is within the map validity range.
-
Int_t getMapBin(Int_t iX, Int_t iY, Int_t iZ)
Find the vector entry of the field map data given the bins iX, iY and iZ.
+
Float_t BInterCalc(CoordAxis theAxis)
+
binPair getBinInfo(Float_t x, CoordAxis theAxis)
Get the bin number and fractional distance from the leftmost bin edge.
+
Bool_t insideRange(Float_t x, Float_t y, Float_t z)
Check to see if a point is within the map validity range.
+
Int_t getMapBin(Int_t iX, Int_t iY, Int_t iZ)
Find the vector entry of the field map data given the bins iX, iY and iZ.
@@ -1186,44 +1189,44 @@

Returns
the bin number and fractional distance from the leftmost bin edge as a pair
-

Definition at line 430 of file ShipBFieldMap.cxx.

-
431 {
-
432 
-
433  Float_t du(0.0), uMin(0.0), Nu(0);
-
434 
-
435  if (theAxis == ShipBFieldMap::xAxis) {
-
436  du = dx_; uMin = xMin_; Nu = Nx_;
-
437  } else if (theAxis == ShipBFieldMap::yAxis) {
-
438  du = dy_; uMin = yMin_; Nu = Ny_;
-
439  } else if (theAxis == ShipBFieldMap::zAxis) {
-
440  du = dz_; uMin = zMin_; Nu = Nz_;
-
441  }
-
442 
-
443  Int_t iBin(-1);
-
444  Float_t fracL(0.0);
+

Definition at line 433 of file ShipBFieldMap.cxx.

+
434 {
+
435 
+
436  Float_t du(0.0), uMin(0.0), Nu(0);
+
437 
+
438  if (theAxis == ShipBFieldMap::xAxis) {
+
439  du = dx_; uMin = xMin_; Nu = Nx_;
+
440  } else if (theAxis == ShipBFieldMap::yAxis) {
+
441  du = dy_; uMin = yMin_; Nu = Ny_;
+
442  } else if (theAxis == ShipBFieldMap::zAxis) {
+
443  du = dz_; uMin = zMin_; Nu = Nz_;
+
444  }
445 
-
446  if (du > 1e-10) {
-
447 
-
448  // Get the number of fractional bin widths the point is from the first volume bin
-
449  Float_t dist = (u - uMin)/du;
-
450  // Get the integer equivalent of this distance, which is the bin number
-
451  iBin = static_cast<Int_t>(dist);
-
452  // Get the actual fractional distance of the point from the leftmost bin edge
-
453  fracL = (dist - iBin*1.0);
-
454 
-
455  }
-
456 
-
457  // Check that the bin number is valid
-
458  if (iBin < 0 || iBin >= Nu) {
-
459  iBin = -1; fracL = 0.0;
-
460  }
-
461 
-
462  // Return the bin number and fractional distance of the point from the leftmost bin edge
-
463  binPair binInfo(iBin, fracL);
+
446  Int_t iBin(-1);
+
447  Float_t fracL(0.0);
+
448 
+
449  if (du > 1e-10) {
+
450 
+
451  // Get the number of fractional bin widths the point is from the first volume bin
+
452  Float_t dist = (u - uMin)/du;
+
453  // Get the integer equivalent of this distance, which is the bin number
+
454  iBin = static_cast<Int_t>(dist);
+
455  // Get the actual fractional distance of the point from the leftmost bin edge
+
456  fracL = (dist - iBin*1.0);
+
457 
+
458  }
+
459 
+
460  // Check that the bin number is valid
+
461  if (iBin < 0 || iBin >= Nu) {
+
462  iBin = -1; fracL = 0.0;
+
463  }
464 
-
465  return binInfo;
-
466 
-
467 }
+
465  // Return the bin number and fractional distance of the point from the leftmost bin edge
+
466  binPair binInfo(iBin, fracL);
+
467 
+
468  return binInfo;
+
469 
+
470 }

@@ -1402,22 +1405,22 @@

Returns
the index entry for the field map data vector
-

Definition at line 469 of file ShipBFieldMap.cxx.

-
470 {
-
471 
-
472  // Get the index of the map entry corresponding to the x,y,z bins.
-
473  // Remember that the map is ordered in ascending z, y, then x
+

Definition at line 472 of file ShipBFieldMap.cxx.

+
473 {
474 
-
475  Int_t index = (iX*Ny_ + iY)*Nz_ + iZ;
-
476  if (index < 0) {
-
477  index = 0;
-
478  } else if (index >= N_) {
-
479  index = N_-1;
-
480  }
-
481 
-
482  return index;
-
483 
-
484 }
+
475  // Get the index of the map entry corresponding to the x,y,z bins.
+
476  // Remember that the map is ordered in ascending z, y, then x
+
477 
+
478  Int_t index = (iX*Ny_ + iY)*Nz_ + iZ;
+
479  if (index < 0) {
+
480  index = 0;
+
481  } else if (index >= N_) {
+
482  index = N_-1;
+
483  }
+
484 
+
485  return index;
+
486 
+
487 }
@@ -2128,34 +2131,34 @@

Definition at line 200 of file ShipBFieldMap.cxx.

-
201 {
-
202 
-
203  if (initialised_ == kFALSE) {
-
204 
-
205  if (isCopy_ == kFALSE) {this->readMapFile();}
-
206 
-
207  // Set the global co-ordinate translation and rotation info
-
208  if (fabs(phi_) > 1e-6 || fabs(theta_) > 1e-6 || fabs(psi_) > 1e-6) {
+

Definition at line 203 of file ShipBFieldMap.cxx.

+
204 {
+
205 
+
206  if (initialised_ == kFALSE) {
+
207 
+
208  if (isCopy_ == kFALSE) {this->readMapFile();}
209 
-
210  // We have non-zero rotation angles. Create a combined translation and rotation
-
211  TGeoTranslation tr("offsets", xOffset_, yOffset_, zOffset_);
-
212  TGeoRotation rot("angles", phi_, theta_, psi_);
-
213  theTrans_ = new TGeoCombiTrans(tr, rot);
-
214 
-
215  } else {
-
216 
-
217  // We only need a translation
-
218  theTrans_ = new TGeoTranslation("offsets", xOffset_, yOffset_, zOffset_);
+
210  // Set the global co-ordinate translation and rotation info
+
211  if (fabs(phi_) > 1e-6 || fabs(theta_) > 1e-6 || fabs(psi_) > 1e-6) {
+
212 
+
213  // We have non-zero rotation angles. Create a combined translation and rotation
+
214  TGeoTranslation tr("offsets", xOffset_, yOffset_, zOffset_);
+
215  TGeoRotation rot("angles", phi_, theta_, psi_);
+
216  theTrans_ = new TGeoCombiTrans(tr, rot);
+
217 
+
218  } else {
219 
-
220  }
-
221 
-
222  initialised_ = kTRUE;
-
223 
-
224  }
-
225 
-
226 }
-
void readMapFile()
Read the field map data and store the information internally.
+
220  // We only need a translation
+
221  theTrans_ = new TGeoTranslation("offsets", xOffset_, yOffset_, zOffset_);
+
222 
+
223  }
+
224 
+
225  initialised_ = kTRUE;
+
226 
+
227  }
+
228 
+
229 }
+
void readMapFile()
Read the field map data and store the information internally.
@@ -2210,17 +2213,17 @@

Returns
true/false if the point is inside the field map range
-

Definition at line 381 of file ShipBFieldMap.cxx.

-
382 {
-
383 
-
384  Bool_t inside(kFALSE);
-
385 
-
386  if (x >= xMin_ && x <= xMax_ && y >= yMin_ &&
-
387  y <= yMax_ && z >= zMin_ && z <= zMax_) {inside = kTRUE;}
+

Definition at line 384 of file ShipBFieldMap.cxx.

+
385 {
+
386 
+
387  Bool_t inside(kFALSE);
388 
-
389  return inside;
-
390 
-
391 }
+
389  if (x >= xMin_ && x <= xMax_ && y >= yMin_ &&
+
390  y <= yMax_ && z >= zMin_ && z <= zMax_) {inside = kTRUE;}
+
391 
+
392  return inside;
+
393 
+
394 }
@@ -2308,26 +2311,26 @@

Definition at line 228 of file ShipBFieldMap.cxx.

-
229 {
-
230 
-
231  std::cout<<"ShipBFieldMap::readMapFile() creating field "<<this->GetName()
-
232  <<" using file "<<mapFileName_<<std::endl;
+

Definition at line 231 of file ShipBFieldMap.cxx.

+
232 {
233 
-
234  // Check to see if we have a ROOT file
-
235  if (mapFileName_.find(".root") != std::string::npos) {
+
234  std::cout<<"ShipBFieldMap::readMapFile() creating field "<<this->GetName()
+
235  <<" using file "<<mapFileName_<<std::endl;
236 
-
237  this->readRootFile();
-
238 
-
239  } else {
-
240 
-
241  this->readTextFile();
-
242 
-
243  }
-
244 
-
245 }
-
void readTextFile()
Process the text file containing the field map data.
-
void readRootFile()
Process the ROOT file containing the field map data.
+
237  // Check to see if we have a ROOT file
+
238  if (mapFileName_.find(".root") != std::string::npos) {
+
239 
+
240  this->readRootFile();
+
241 
+
242  } else {
+
243 
+
244  this->readTextFile();
+
245 
+
246  }
+
247 
+
248 }
+
void readTextFile()
Process the text file containing the field map data.
+
void readRootFile()
Process the ROOT file containing the field map data.
@@ -2356,92 +2359,92 @@

Definition at line 247 of file ShipBFieldMap.cxx.

-
247  {
-
248 
-
249  TFile* theFile = TFile::Open(mapFileName_.c_str());
-
250 
-
251  if (!theFile) {
-
252  std::cout<<"ShipBFieldMap: could not find the file "<<mapFileName_<<std::endl;
-
253  return;
-
254  }
-
255 
-
256  // Coordinate ranges
-
257  TTree* rTree = dynamic_cast<TTree*>(theFile->Get("Range"));
-
258  if (!rTree) {
-
259  std::cout<<"ShipBFieldMap: could not find Range tree in "<<mapFileName_<<std::endl;
-
260  return;
-
261  }
-
262 
-
263  rTree->SetBranchAddress("xMin", &xMin_);
-
264  rTree->SetBranchAddress("xMax", &xMax_);
-
265  rTree->SetBranchAddress("dx", &dx_);
-
266  rTree->SetBranchAddress("yMin", &yMin_);
-
267  rTree->SetBranchAddress("yMax", &yMax_);
-
268  rTree->SetBranchAddress("dy", &dy_);
-
269  rTree->SetBranchAddress("zMin", &zMin_);
-
270  rTree->SetBranchAddress("zMax", &zMax_);
-
271  rTree->SetBranchAddress("dz", &dz_);
-
272 
-
273  // Fill the ranges
-
274  rTree->GetEntry(0);
+

Definition at line 250 of file ShipBFieldMap.cxx.

+
250  {
+
251 
+
252  TFile* theFile = TFile::Open(mapFileName_.c_str());
+
253 
+
254  if (!theFile) {
+
255  std::cout<<"ShipBFieldMap: could not find the file "<<mapFileName_<<std::endl;
+
256  return;
+
257  }
+
258 
+
259  // Coordinate ranges
+
260  TTree* rTree = dynamic_cast<TTree*>(theFile->Get("Range"));
+
261  if (!rTree) {
+
262  std::cout<<"ShipBFieldMap: could not find Range tree in "<<mapFileName_<<std::endl;
+
263  return;
+
264  }
+
265 
+
266  rTree->SetBranchAddress("xMin", &xMin_);
+
267  rTree->SetBranchAddress("xMax", &xMax_);
+
268  rTree->SetBranchAddress("dx", &dx_);
+
269  rTree->SetBranchAddress("yMin", &yMin_);
+
270  rTree->SetBranchAddress("yMax", &yMax_);
+
271  rTree->SetBranchAddress("dy", &dy_);
+
272  rTree->SetBranchAddress("zMin", &zMin_);
+
273  rTree->SetBranchAddress("zMax", &zMax_);
+
274  rTree->SetBranchAddress("dz", &dz_);
275 
-
276  this->setLimits();
-
277 
-
278  // Make sure we don't have a copy
-
279  if (isCopy_ == kFALSE) {
+
276  // Fill the ranges
+
277  rTree->GetEntry(0);
+
278 
+
279  this->setLimits();
280 
-
281  // The data is expected to contain Bx,By,Bz data values
-
282  // in ascending z,y,x co-ordinate order
+
281  // Make sure we don't have a copy
+
282  if (isCopy_ == kFALSE) {
283 
-
284  TTree* dTree = dynamic_cast<TTree*>(theFile->Get("Data"));
-
285  if (!dTree) {
-
286  std::cout<<"ShipBFieldMap: could not find Data tree in "<<mapFileName_<<std::endl;
-
287  return;
-
288  }
-
289 
-
290  Float_t Bx, By, Bz;
-
291  // Only enable the field components
-
292  dTree->SetBranchStatus("*", 0);
-
293  dTree->SetBranchStatus("Bx", 1);
-
294  dTree->SetBranchStatus("By", 1);
-
295  dTree->SetBranchStatus("Bz", 1);
-
296 
-
297  dTree->SetBranchAddress("Bx", &Bx);
-
298  dTree->SetBranchAddress("By", &By);
-
299  dTree->SetBranchAddress("Bz", &Bz);
-
300 
-
301  Int_t nEntries = dTree->GetEntries();
-
302  if (nEntries != N_) {
-
303  std::cout<<"Expected "<<N_<<" field map entries but found "<<nEntries<<std::endl;
-
304  nEntries = 0;
-
305  }
-
306 
-
307  fieldMap_->reserve(nEntries);
-
308 
-
309  for (Int_t i = 0; i < nEntries; i++) {
-
310 
-
311  dTree->GetEntry(i);
-
312 
-
313  // B field values are in Tesla. This means these values are multiplied
-
314  // by a factor of ten since both FairRoot and the VMC interface use kGauss
-
315  Bx *= Tesla_;
-
316  By *= Tesla_;
-
317  Bz *= Tesla_;
-
318 
-
319  // Store the B field 3-vector
-
320  std::vector<Float_t> BVector(3);
-
321  BVector[0] = Bx; BVector[1] = By; BVector[2] = Bz;
-
322  fieldMap_->push_back(BVector);
-
323 
-
324  }
-
325 
-
326  }
-
327 
-
328  theFile->Close();
-
329 
-
330 }
- +
284  // The data is expected to contain Bx,By,Bz data values
+
285  // in ascending z,y,x co-ordinate order
+
286 
+
287  TTree* dTree = dynamic_cast<TTree*>(theFile->Get("Data"));
+
288  if (!dTree) {
+
289  std::cout<<"ShipBFieldMap: could not find Data tree in "<<mapFileName_<<std::endl;
+
290  return;
+
291  }
+
292 
+
293  Float_t Bx, By, Bz;
+
294  // Only enable the field components
+
295  dTree->SetBranchStatus("*", 0);
+
296  dTree->SetBranchStatus("Bx", 1);
+
297  dTree->SetBranchStatus("By", 1);
+
298  dTree->SetBranchStatus("Bz", 1);
+
299 
+
300  dTree->SetBranchAddress("Bx", &Bx);
+
301  dTree->SetBranchAddress("By", &By);
+
302  dTree->SetBranchAddress("Bz", &Bz);
+
303 
+
304  Int_t nEntries = dTree->GetEntries();
+
305  if (nEntries != N_) {
+
306  std::cout<<"Expected "<<N_<<" field map entries but found "<<nEntries<<std::endl;
+
307  nEntries = 0;
+
308  }
+
309 
+
310  fieldMap_->reserve(nEntries);
+
311 
+
312  for (Int_t i = 0; i < nEntries; i++) {
+
313 
+
314  dTree->GetEntry(i);
+
315 
+
316  // B field values are in Tesla. This means these values are multiplied
+
317  // by a factor of ten since both FairRoot and the VMC interface use kGauss
+
318  Bx *= Tesla_;
+
319  By *= Tesla_;
+
320  Bz *= Tesla_;
+
321 
+
322  // Store the B field 3-vector
+
323  std::vector<Float_t> BVector(3);
+
324  BVector[0] = Bx; BVector[1] = By; BVector[2] = Bz;
+
325  fieldMap_->push_back(BVector);
+
326 
+
327  }
+
328 
+
329  }
+
330 
+
331  theFile->Close();
+
332 
+
333 }
+
int i
Definition: ShipAna.py:86
@@ -2475,55 +2478,55 @@

Definition at line 332 of file ShipBFieldMap.cxx.

-
332  {
-
333 
-
334  std::ifstream getData(mapFileName_.c_str());
-
335 
-
336  std::string tmpString("");
-
337 
-
338  getData >> tmpString >> xMin_ >> xMax_ >> dx_
-
339  >> yMin_ >> yMax_ >> dy_ >> zMin_ >> zMax_ >> dz_;
+

Definition at line 335 of file ShipBFieldMap.cxx.

+
335  {
+
336 
+
337  std::ifstream getData(mapFileName_.c_str());
+
338 
+
339  std::string tmpString("");
340 
-
341  this->setLimits();
-
342 
-
343  // Check to see if this object is a "copy"
-
344  if (isCopy_ == kFALSE) {
+
341  getData >> tmpString >> xMin_ >> xMax_ >> dx_
+
342  >> yMin_ >> yMax_ >> dy_ >> zMin_ >> zMax_ >> dz_;
+
343 
+
344  this->setLimits();
345 
-
346  // Read the field map and store the magnetic field components
-
347 
-
348  // Second line can be ignored since they are
-
349  // just labels for the data columns for readability
-
350  getData >> tmpString >> tmpString >> tmpString;
-
351 
-
352  // The remaining lines contain Bx,By,Bz data values
-
353  // in ascending z,y,x co-ord order
-
354  fieldMap_->reserve(N_);
-
355 
-
356  Float_t Bx(0.0), By(0.0), Bz(0.0);
-
357 
-
358  for (Int_t i = 0; i < N_; i++) {
-
359 
-
360  getData >> Bx >> By >> Bz;
-
361 
-
362  // B field values are in Tesla. This means these values are multiplied
-
363  // by a factor of ten since both FairRoot and the VMC interface use kGauss
-
364  Bx *= Tesla_;
-
365  By *= Tesla_;
-
366  Bz *= Tesla_;
-
367 
-
368  // Store the B field 3-vector
-
369  std::vector<Float_t> BVector(3);
-
370  BVector[0] = Bx; BVector[1] = By; BVector[2] = Bz;
-
371  fieldMap_->push_back(BVector);
-
372 
-
373  }
-
374 
-
375  }
-
376 
-
377  getData.close();
-
378 
-
379 }
+
346  // Check to see if this object is a "copy"
+
347  if (isCopy_ == kFALSE) {
+
348 
+
349  // Read the field map and store the magnetic field components
+
350 
+
351  // Second line can be ignored since they are
+
352  // just labels for the data columns for readability
+
353  getData >> tmpString >> tmpString >> tmpString;
+
354 
+
355  // The remaining lines contain Bx,By,Bz data values
+
356  // in ascending z,y,x co-ord order
+
357  fieldMap_->reserve(N_);
+
358 
+
359  Float_t Bx(0.0), By(0.0), Bz(0.0);
+
360 
+
361  for (Int_t i = 0; i < N_; i++) {
+
362 
+
363  getData >> Bx >> By >> Bz;
+
364 
+
365  // B field values are in Tesla. This means these values are multiplied
+
366  // by a factor of ten since both FairRoot and the VMC interface use kGauss
+
367  Bx *= Tesla_;
+
368  By *= Tesla_;
+
369  Bz *= Tesla_;
+
370 
+
371  // Store the B field 3-vector
+
372  std::vector<Float_t> BVector(3);
+
373  BVector[0] = Bx; BVector[1] = By; BVector[2] = Bz;
+
374  fieldMap_->push_back(BVector);
+
375 
+
376  }
+
377 
+
378  }
+
379 
+
380  getData.close();
+
381 
+
382 }

@@ -2550,42 +2553,42 @@

-

Definition at line 393 of file ShipBFieldMap.cxx.

-
393  {
-
394 
-
395  // Since the default SHIP distance unit is cm, we do not need to convert
-
396  // these map limits, i.e. cm = 1 already
+

Definition at line 396 of file ShipBFieldMap.cxx.

+
396  {
397 
-
398  xRange_ = xMax_ - xMin_;
-
399  yRange_ = yMax_ - yMin_;
-
400  zRange_ = zMax_ - zMin_;
-
401 
-
402  // Find the number of bins using the limits and bin sizes. The number of bins
-
403  // includes both the minimum and maximum values. To ensure correct rounding
-
404  // up to the nearest integer we need to add 1.5 not 1.0.
-
405  if (dx_ > 0.0) {
-
406  Nx_ = static_cast<Int_t>(((xMax_ - xMin_)/dx_) + 1.5);
-
407  }
-
408  if (dy_ > 0.0) {
-
409  Ny_ = static_cast<Int_t>(((yMax_ - yMin_)/dy_) + 1.5);
+
398  // Since the default SHIP distance unit is cm, we do not need to convert
+
399  // these map limits, i.e. cm = 1 already
+
400 
+
401  xRange_ = xMax_ - xMin_;
+
402  yRange_ = yMax_ - yMin_;
+
403  zRange_ = zMax_ - zMin_;
+
404 
+
405  // Find the number of bins using the limits and bin sizes. The number of bins
+
406  // includes both the minimum and maximum values. To ensure correct rounding
+
407  // up to the nearest integer we need to add 1.5 not 1.0.
+
408  if (dx_ > 0.0) {
+
409  Nx_ = static_cast<Int_t>(((xMax_ - xMin_)/dx_) + 1.5);
410  }
-
411  if (dz_ > 0.0) {
-
412  Nz_ = static_cast<Int_t>(((zMax_ - zMin_)/dz_) + 1.5);
+
411  if (dy_ > 0.0) {
+
412  Ny_ = static_cast<Int_t>(((yMax_ - yMin_)/dy_) + 1.5);
413  }
-
414 
-
415  N_ = Nx_*Ny_*Nz_;
-
416 
-
417  std::cout<<"x limits: "<<xMin_<<", "<<xMax_<<", dx = "<<dx_<<std::endl;
-
418  std::cout<<"y limits: "<<yMin_<<", "<<yMax_<<", dy = "<<dy_<<std::endl;
-
419  std::cout<<"z limits: "<<zMin_<<", "<<zMax_<<", dz = "<<dz_<<std::endl;
-
420 
-
421  std::cout<<"Offsets: x = "<<xOffset_<<", y = "<<yOffset_<<", z = "<<zOffset_<<std::endl;
-
422  std::cout<<"Angles : phi = "<<phi_<<", theta = "<<theta_<<", psi = "<<psi_<<std::endl;
+
414  if (dz_ > 0.0) {
+
415  Nz_ = static_cast<Int_t>(((zMax_ - zMin_)/dz_) + 1.5);
+
416  }
+
417 
+
418  N_ = Nx_*Ny_*Nz_;
+
419 
+
420  std::cout<<"x limits: "<<xMin_<<", "<<xMax_<<", dx = "<<dx_<<std::endl;
+
421  std::cout<<"y limits: "<<yMin_<<", "<<yMax_<<", dy = "<<dy_<<std::endl;
+
422  std::cout<<"z limits: "<<zMin_<<", "<<zMax_<<", dz = "<<dz_<<std::endl;
423 
-
424  std::cout<<"Total number of bins = "<<N_
-
425  <<"; Nx = "<<Nx_<<", Ny = "<<Ny_<<", Nz = "<<Nz_<<std::endl;
+
424  std::cout<<"Offsets: x = "<<xOffset_<<", y = "<<yOffset_<<", z = "<<zOffset_<<std::endl;
+
425  std::cout<<"Angles : phi = "<<phi_<<", theta = "<<theta_<<", psi = "<<psi_<<std::endl;
426 
-
427 }
+
427  std::cout<<"Total number of bins = "<<N_
+
428  <<"; Nx = "<<Nx_<<", Ny = "<<Ny_<<", Nz = "<<Nz_<<std::endl;
+
429 
+
430 }

diff --git a/md_CHANGELOG.html b/md_CHANGELOG.html index 8a57be297..ed5bffe0d 100644 --- a/md_CHANGELOG.html +++ b/md_CHANGELOG.html @@ -96,6 +96,7 @@

Fixed

  • Use ConstructedAt instead of remove pythonization for TClonesArray
  • +
  • Octant symmetry was incorrect for B_z when using field maps (reported and fixed by M. Ferro-Luzzi)

Changed