STEER/AliMagFC.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 //     Constant magnetic field class
  20 //     Used by AliRun class
  21 //     Author:
  22 //-------------------------------------------------------------------------
  23
  24 #include <stdlib.h>
  25
  26 #include "AliLog.h"
  27 #include "AliMagFC.h"
  28
  29 ClassImp(AliMagFC)
  30
  31 //________________________________________
  32 AliMagFC::AliMagFC()
  33     :AliMagF(),
  34     fCompensator(kFALSE),
  35     fBeamType(kBeamTypepp),
  36     fBeamEnergy(0),
  37     fQuadGradient(0),
  38     fDipoleField(0),
  39     fCCorrField(0),
  40     fACorr1Field(0),
  41     fACorr2Field(0)
  42 {
  43   //
  44   // Default constructor
  45   //
  46 }
  47
  48 //________________________________________
  49 AliMagFC::AliMagFC(const char *name, const char *title, Int_t integ,
  50                    Float_t factor, Float_t fmax)
  51     : AliMagF(name,title,integ,factor,fmax),
  52     fCompensator(kFALSE),
  53     fBeamType(kBeamTypepp),
  54     fBeamEnergy(7000.),
  55     fQuadGradient(0),
  56     fDipoleField(0),
  57     fCCorrField(0),
  58     fACorr1Field(0),
  59     fACorr2Field(0)
  60
  61 {
  62   //
  63   // Standard constructor
  64   //
  65   fType = kConst;
  66   fMap  = 1;
  67
  68 }
  69
  70 //________________________________________
  71 void AliMagFC::Field(Float_t *x, Float_t *b) const
  72 {
  73   //
  74   // Method to return the field in a point
  75   //
  76   b[0]=b[1]=b[2]=0;
  77   if(fMap==1) {
  78     if(TMath::Abs(x[2])<700 && x[0]*x[0]+(x[1]+30)*(x[1]+30) < 560*560) {
  79       b[2]=2;
  80     }
  81     else {
  82       if(-725 >= x[2] && x[2] >= -1225 ){
  83         Float_t dz = TMath::Abs(-975-x[2])*0.01;
  84         b[0] = - (1-0.1*dz*dz)*7;
  85         if(fFactor!=1) {
  86             b[0]*=fFactor;
  87             b[1]*=fFactor;
  88             b[2]*=fFactor;
  89         }
  90       }
  91       else {
  92           ZDCField(x, b);
  93       }
  94     }
  95
  96   }
  97   else {
  98       AliFatal(Form("Invalid field map for constant field %d",fMap));
  99   }
 100 }
 101
 102 //___________________________________________________
 103 void AliMagFC::ZDCField(Float_t *x, Float_t *b) const
 104 {
 105   // ---- This is the ZDC part
 106
 107   Float_t rad2 = x[0] * x[0] + x[1] * x[1];
 108   static Bool_t init = kFALSE;
 109
 110   if (! init) {
 111       init = kTRUE;
 112       //////////////////////////////////////////////////////////////////////
 113       // ---- Magnetic field values (according to beam type and energy) ----
 114       if(fBeamType==kBeamTypepp && fBeamEnergy == 5000.){
 115           // p-p @ 5+5 TeV
 116           fQuadGradient = 15.7145;
 117           fDipoleField  = 27.0558;
 118           // SIDE C
 119           fCCorrField   = 9.7017;
 120           // SIDE A
 121           fACorr1Field  = -13.2143;
 122           fACorr2Field  = -11.9909;
 123       } else if (fBeamType == kBeamTypepp && fBeamEnergy == 450.) {
 124           // p-p 0.45+0.45 TeV
 125           Float_t const kEnergyRatio = fBeamEnergy / 7000.;
 126
 127           fQuadGradient = 22.0002 * kEnergyRatio;
 128           fDipoleField  = 37.8781 * kEnergyRatio;
 129           // SIDE C
 130           fCCorrField   =  9.6908;
 131           // SIDE A
 132           fACorr1Field  = -13.2014;
 133           fACorr2Field  = -9.6908;
 134       } else if ((fBeamType == kBeamTypepp && fBeamEnergy == 7000.) ||
 135                  (fBeamType == kBeamTypeAA))
 136       {
 137           // Pb-Pb @ 2.7+2.7 TeV or p-p @ 7+7 TeV
 138           fQuadGradient = 22.0002;
 139           fDipoleField  = 37.8781;
 140           // SIDE C
 141           fCCorrField   = 9.6908;
 142           // SIDE A
 143           fACorr1Field  = -13.2014;
 144           fACorr2Field  = -9.6908;
 145       }
 146
 147       printf("Machine field %5d %13.3f %13.3f \n", fBeamType, fBeamEnergy, fDipoleField);
 148   }
 149
 150
 151   // SIDE C **************************************************
 152   if(x[2]<0.){
 153     if(x[2] < kCCorrBegin && x[2] > kCCorrEnd && rad2 < kCCorrSqRadius){
 154         if (fFactor != 0.) {
 155             b[0] = fCCorrField;
 156             b[1] = 0.;
 157             b[2] = 0.;
 158         }
 159     }
 160     else if(x[2] < kCQ1Begin && x[2] > kCQ1End && rad2 < kCQ1SqRadius){
 161         b[0] = fQuadGradient*x[1];
 162         b[1] = fQuadGradient*x[0];
 163         b[2] = 0.;
 164     }
 165     else if(x[2] < kCQ2Begin && x[2] > kCQ2End && rad2 < kCQ2SqRadius){
 166         b[0] = -fQuadGradient*x[1];
 167         b[1] = -fQuadGradient*x[0];
 168         b[2] = 0.;
 169     }
 170     else if(x[2] < kCQ3Begin && x[2] > kCQ3End && rad2 < kCQ3SqRadius){
 171         b[0] = -fQuadGradient*x[1];
 172         b[1] = -fQuadGradient*x[0];
 173         b[2] = 0.;
 174     }
 175     else if(x[2] < kCQ4Begin && x[2] > kCQ4End && rad2 < kCQ4SqRadius){
 176         b[0] = fQuadGradient*x[1];
 177         b[1] = fQuadGradient*x[0];
 178         b[2] = 0.;
 179     }
 180     else if(x[2] < kCD1Begin && x[2] > kCD1End && rad2 < kCD1SqRadius){
 181         b[1] = fDipoleField;
 182         b[2] = 0.;
 183         b[2] = 0.;
 184     }
 185     else if(x[2] < kCD2Begin && x[2] > kCD2End){
 186         if(((x[0]-kCD2XCentre1)*(x[0]-kCD2XCentre1)+(x[1]*x[1]))<kCD2SqRadius
 187            || ((x[0]-kCD2XCentre2)*(x[0]-kCD2XCentre2)+(x[1]*x[1]))<kCD2SqRadius){
 188           b[1] = -fDipoleField;
 189           b[2] = 0.;
 190           b[2] = 0.;
 191         }
 192     }
 193   }
 194
 195   // SIDE A **************************************************
 196   else{
 197     if(fCompensator && (x[2] > kACorr1Begin && x[2] < kACorr1End) && rad2 < kCCorr1SqRadius) {
 198       // Compensator magnet at z = 1075 m
 199         if (fFactor != 0.) {
 200             b[0] = fACorr1Field;
 201             b[1] = 0.;
 202             b[2] = 0.;
 203         }
 204         return;
 205     }
 206
 207     if(x[2] > kACorr2Begin && x[2] < kACorr2End && rad2 < kCCorr2SqRadius){
 208         if (fFactor != 0.) {
 209             b[0] = fACorr2Field;
 210             b[1] = 0.;
 211             b[2] = 0.;
 212         }
 213     }
 214     else if(x[2] > kAQ1Begin && x[2] < kAQ1End && rad2 < kAQ1SqRadius){
 215         // First quadrupole of inner triplet de-focussing in x-direction
 216         b[0] = -fQuadGradient*x[1];
 217         b[1] = -fQuadGradient*x[0];
 218         b[2] = 0.;
 219     }
 220     else if(x[2] > kAQ2Begin && x[2] < kAQ2End && rad2 < kAQ2SqRadius){
 221         b[0] = fQuadGradient*x[1];
 222         b[1] = fQuadGradient*x[0];
 223         b[2] = 0.;
 224     }
 225     else if(x[2] > kAQ3Begin && x[2] < kAQ3End && rad2 < kAQ3SqRadius){
 226         b[0] = fQuadGradient*x[1];
 227         b[1] = fQuadGradient*x[0];
 228         b[2] = 0.;
 229     }
 230     else if(x[2] > kAQ4Begin && x[2] < kAQ4End && rad2 < kAQ4SqRadius){
 231         b[0] = -fQuadGradient*x[1];
 232         b[1] = -fQuadGradient*x[0];
 233         b[2] = 0.;
 234     }
 235     else if(x[2] > kAD1Begin && x[2] < kAD1End && rad2 < kAD1SqRadius){
 236         b[0] = 0.;
 237         b[1] = -fDipoleField;
 238         b[2] = 0.;
 239     }
 240     else if(x[2] > kAD2Begin && x[2] < kAD2End){
 241         if(((x[0]-kAD2XCentre1)*(x[0]-kAD2XCentre1)+(x[1]*x[1])) < kAD2SqRadius
 242            || ((x[0]-kAD2XCentre2)*(x[0]-kAD2XCentre2)+(x[1]*x[1])) < kAD2SqRadius){
 243             b[1] = fDipoleField;
 244         }
 245     }
 246   }
 247 }
 248