BLCMDhelicaldipole Class Reference

Inheritance diagram for BLCMDhelicaldipole:

BLElement BLManager::RunAction BLCommand

List of all members.


Detailed Description

BLCMDhelicaldipole implements helical dipole magnet with a cylindrical field volume.

model 1 The magnetic field is an ideal dipole that rotates helically along z. This is a pure field -- there is no iron or physical volume. model 2 The field is determined from solutions of cylindrical Maxwellian. RPJ and YSD designed this. (MuNote0284) model 3 The field is determined from solutions of cylindrical Maxwellian. This time you take into account the boundary condition. T. Tominaka et al designed this. (NIM A459:398) model 4 The field strength can be modulated as z.

Public Member Functions

 BLCMDhelicaldipole ()
 Default constructor. Defines the command, args, etc.
virtual ~BLCMDhelicaldipole ()
 Destructor.
 BLCMDhelicaldipole (const BLCMDhelicaldipole &r)
 Copy constructor.
BLElementclone ()
 clone()
G4String commandName ()
 commandName() returns "helicaldipole".
int command (BLArgumentVector &argv, BLArgumentMap &namedArgs)
 command() implements the helicaldipole command.
void defineNamedArgs ()
 defineNamedArgs() defines the named arguments for the command.
void construct (G4RotationMatrix *relativeRotation, G4ThreeVector relativePosition, G4LogicalVolume *parent, G4String parentName, G4RotationMatrix *parentRotation, G4ThreeVector parentPosition)
 construct() - construct the helicaldipole magnet
G4double getLength ()
 getLength() returns the fieldLength of the hd
G4double getWidth ()
 getWidth() returns the outer radius of the hd
G4double getHeight ()
 getHeight() returns the outer radius of the hd
G4bool isOK ()
 isOK() returns true.
bool isOutside (G4ThreeVector &local, G4double tolerance)
 isOutside() from BLElement.
void generatePoints (int npoints, std::vector< G4ThreeVector > &v)
 generatePoints() from BLElement.
void BeginOfRunAction (const G4Run *run)
 BeginOfRunAction() from BLManager::RunAction.
void EndOfRunAction (const G4Run *run)
 EndOfRunAction() from BLManager::RunAction.

Private Attributes

G4double radius
G4double length
G4double b
G4double bprime
G4double bpp
G4double rr
G4double lambda
G4double model
G4double psi0
G4double phi0
G4double Bsolenoid
G4double ez
BLCoordinateTransform global2local
G4Polymarker markers

Friends

class HelicalDipoleField


Constructor & Destructor Documentation

BLCMDhelicaldipole::BLCMDhelicaldipole (  ) 

Default constructor. Defines the command, args, etc.

References b, BLCMDTYPE_ELEMENT, bpp, bprime, Bsolenoid, ez, lambda, length, model, phi0, psi0, radius, BLCommand::registerCommand(), rr, BLCommand::setDescription(), and BLCommand::setSynopsis().

Referenced by clone(), and command().

00197                                        : BLElement(), BLManager::RunAction()
00198 {
00199         // register the commandName(), and its synopsis and description.
00200         registerCommand(BLCMDTYPE_ELEMENT);
00201         setSynopsis("construct a helicaldipole magnet.");
00202         setDescription(
00203                 "The field region is a cylinder with a helical dipole\n"
00204                 "field plus a solenoidial field.  The simple model=1 \n"
00205                 "provides just a sine and cosine transverse dependence,\n"
00206                 "while the maxwellian model=2 has both dipole and quadrupole\n"
00207                 "terms.  Both the dipole scale bD [T] and quadrupole scale bQ [T/m]\n"
00208                 "are now at rho=0; the user must determine the correct values externally.\n\n"
00209                 "Note that this Element generates a magnetic field only,\n"
00210                 "and only within the cylinder defined by length and radius.\n"
00211                 "So it has no solid associated with it, and is invisible.\n");
00212 
00213         // provide initial values for fields
00214         radius = 0.0;
00215         length = 0.0;
00216         b = 0.0;
00217         bprime = 0.0;
00218         bpp    = 0.0;
00219         rr     = 159.155;
00220         model = 3. ;
00221         lambda = 0.0;
00222         phi0 = 0.0;
00223         psi0 = 0.0;
00224         Bsolenoid = 0.0;
00225         ez = 0.0;
00226 }

virtual BLCMDhelicaldipole::~BLCMDhelicaldipole (  )  [inline, virtual]

Destructor.

00113 { }

BLCMDhelicaldipole::BLCMDhelicaldipole ( const BLCMDhelicaldipole r  ) 

Copy constructor.

References b, bpp, bprime, Bsolenoid, ez, lambda, length, model, phi0, psi0, radius, and rr.

00229                                                                   : BLElement(r), 
00230                                                         BLManager::RunAction(r)
00231 {
00232         // copy fields one at a time (transfers default values from the
00233         // default object to this new object).
00234         radius = r.radius;
00235         length = r.length;
00236         b = r.b;
00237         bprime = r.bprime;
00238         bpp    = r.bpp;
00239         rr     = r.rr;
00240         lambda = r.lambda;
00241         phi0 = r.phi0;
00242         model = r.model;
00243         Bsolenoid = r.Bsolenoid;
00244         ez = r.ez;
00245         psi0 = r.psi0;
00246 }


Member Function Documentation

BLElement* BLCMDhelicaldipole::clone (  )  [inline, virtual]

clone()

Implements BLElement.

References BLCMDhelicaldipole().

00119 { return new BLCMDhelicaldipole(*this); }

G4String BLCMDhelicaldipole::commandName (  )  [inline, virtual]

commandName() returns "helicaldipole".

Implements BLCommand.

00122 { return "helicaldipole"; }

int BLCMDhelicaldipole::command ( BLArgumentVector argv,
BLArgumentMap namedArgs 
) [virtual]

command() implements the helicaldipole command.

Implements BLCommand.

References BLCMDhelicaldipole(), BLCommand::handleNamedArgs(), BLCommand::print(), BLCommand::printError(), and BLElement::setName().

00249 {
00250         const char *CoDE="BLCMDhelicaldipole::command";
00251 
00252         if(argv.size() != 1) {
00253                 printError("helicaldipole: Invalid command, must have name");
00254                 return -1;
00255         }
00256 
00257         if(argv[0] == "default") {
00258                 return defaultHelicalDipole.handleNamedArgs(namedArgs);
00259         }
00260 
00261         BLCMDhelicaldipole *t = new BLCMDhelicaldipole(defaultHelicalDipole);
00262         t->setName(argv[0]);
00263         int retval = t->handleNamedArgs(namedArgs);
00264 
00265         t->print(argv[0]);
00266 
00267         return retval;
00268 }

void BLCMDhelicaldipole::defineNamedArgs (  )  [virtual]

defineNamedArgs() defines the named arguments for the command.

Reimplemented from BLCommand.

References BLCommand::argDouble(), b, bpp, bprime, Bsolenoid, ez, lambda, length, model, phi0, psi0, radius, and rr.

00271 {
00272         argDouble(radius,"radius","The radius of the field region (mm)",mm);
00273         argDouble(model,"model","The model of field calculated(simple=1, rpj and ysd model=2, ttominaka et al model=3), modulations in bd,bq,bz= 4");
00274         argDouble(length,"length","The length of the field region (mm)",mm);
00275         argDouble(b,"bD","The dipole magnitude at rho=0 (Tesla).",tesla);
00276         argDouble(lambda,"lambda","Helix period along the Z axis (mm).",mm);
00277         argDouble(phi0,"phi0","The phase of the XY field at the entrance (deg).",deg);
00278         argDouble(Bsolenoid,"Bsolenoid","The value of Bsolenoid (Tesla).",tesla);
00279         argDouble(bprime,"bQ","The quadrupole magnitude at rho=0 (Tesla).",tesla/meter);
00280         argDouble(bpp,"bs","The sextupole magnitude at rho=0 (Tesla).",tesla/meter/meter);
00281         argDouble(rr,"rr","Reference radius (mm)",mm);
00282         argDouble(psi0,"psi0","The offset between the dipole term and the quadrupole term (Degrees).",deg);
00283         argDouble(ez,"ez","The base electric field inside the helix channel (GV/m).",1000*megavolt/meter);
00284 }

void BLCMDhelicaldipole::construct ( G4RotationMatrix *  relativeRotation,
G4ThreeVector  relativePosition,
G4LogicalVolume *  parent,
G4String  parentName,
G4RotationMatrix *  parentRotation,
G4ThreeVector  parentPosition 
) [virtual]

construct() - construct the helicaldipole magnet

Implements BLElement.

References BLGlobalField::addElementField(), getLength(), BLElement::getName(), BLManager::getObject(), BLGlobalField::getObject(), global2local, HelicalDipoleField, and BLManager::registerRunAction().

00293 {
00294         G4String thisname = parentName+getName();
00295 
00296         // get globalRotation and globalPosition
00297         G4RotationMatrix *globalRotation = 0;
00298         if(relativeRotation && parentRotation) {
00299                 globalRotation = 
00300                     new G4RotationMatrix(*parentRotation * *relativeRotation);
00301         } else if(relativeRotation) {
00302                 globalRotation = relativeRotation;
00303         } else if(parentRotation) {
00304                 globalRotation = parentRotation;
00305         }
00306         G4ThreeVector globalPosition(relativePosition + parentPosition);
00307         if(parentRotation)
00308                 globalPosition = *parentRotation * relativePosition +
00309                                 parentPosition;
00310 
00311         global2local = BLCoordinateTransform(globalRotation,globalPosition);
00312 
00313         G4double zmin = globalPosition[2]-getLength()/2.0;
00314         G4double zmax = globalPosition[2]+getLength()/2.0;
00315 
00316         HelicalDipoleField *p = new HelicalDipoleField(global2local,this);
00317         BLGlobalField::getObject()->addElementField(p);
00318 
00319         printf("BLCMDhelicaldipole::Construct %s parent=%s relZ=%.1f globZ=%.1f\n"
00320                         "\tzmin=%.1f zmax=%.1f\n",
00321                 thisname.c_str(),parentName.c_str(),relativePosition[2],
00322                 globalPosition[2], zmin,zmax);
00323 
00324         BLManager::getObject()->registerRunAction(this,false);
00325 }

G4double BLCMDhelicaldipole::getLength (  )  [inline, virtual]

getLength() returns the fieldLength of the hd

Implements BLElement.

References length.

Referenced by construct().

00139 { return length; }

G4double BLCMDhelicaldipole::getWidth (  )  [inline, virtual]

getWidth() returns the outer radius of the hd

Implements BLElement.

References radius.

00142 { return radius*2.0; }

G4double BLCMDhelicaldipole::getHeight (  )  [inline, virtual]

getHeight() returns the outer radius of the hd

Implements BLElement.

References radius.

00145 { return radius*2.0; }

G4bool BLCMDhelicaldipole::isOK (  )  [inline, virtual]

isOK() returns true.

Implements BLElement.

00148 { return true; }

bool BLCMDhelicaldipole::isOutside ( G4ThreeVector &  local,
G4double  tolerance 
) [inline, virtual]

isOutside() from BLElement.

Implements BLElement.

00152                 { return true; }

void BLCMDhelicaldipole::generatePoints ( int  npoints,
std::vector< G4ThreeVector > &  v 
) [inline, virtual]

generatePoints() from BLElement.

Implements BLElement.

00156                 { v.clear(); }

void BLCMDhelicaldipole::BeginOfRunAction ( const G4Run *  run  )  [virtual]

BeginOfRunAction() from BLManager::RunAction.

Implements BLManager::RunAction.

References markers.

00328 {
00329         markers.clear();
00330 }

void BLCMDhelicaldipole::EndOfRunAction ( const G4Run *  run  )  [virtual]

EndOfRunAction() from BLManager::RunAction.

Implements BLManager::RunAction.

References BLCoordinateTransform::getGlobal(), global2local, lambda, length, MARKER_SIZE, markers, phi0, and radius.

00333 {
00334         G4VVisManager* pVVisManager = G4VVisManager::GetConcreteInstance();
00335         if (!pVVisManager) return;
00336 
00337 #ifdef STUB // omit markers
00338         double dz = lambda/10.0;
00339         int n = (int)(length/dz) + 1;
00340         for(int i=0; i<n; ++i) {
00341                 G4double local[4], global[4];
00342                 G4double phi = phi0 + i * dz * 2.0*pi/lambda;
00343                 local[0] = radius/2.0*cos(phi);
00344                 local[1] = radius/2.0*sin(phi);
00345                 local[2] = -length/2.0 + i * dz;
00346                 local[3] = 0.0;
00347                 global2local.getGlobal(local,global);
00348                 G4ThreeVector point(global[0],global[1],global[2]);
00349                 markers.push_back(point);
00350         }
00351         markers.SetMarkerType(G4Polymarker::circles);
00352         markers.SetScreenSize(MARKER_SIZE);
00353         markers.SetFillStyle(G4VMarker::filled);
00354         G4VisAttributes va(G4Colour(1.,1.,1.));  // white
00355         markers.SetVisAttributes(&va);
00356 
00357         pVVisManager->Draw(markers);
00358 #endif // STUB
00359 }


Friends And Related Function Documentation

friend class HelicalDipoleField [friend]

Referenced by construct().


Member Data Documentation

G4double BLCMDhelicaldipole::radius [private]

G4double BLCMDhelicaldipole::length [private]

G4double BLCMDhelicaldipole::b [private]

G4double BLCMDhelicaldipole::bprime [private]

G4double BLCMDhelicaldipole::bpp [private]

G4double BLCMDhelicaldipole::rr [private]

G4double BLCMDhelicaldipole::lambda [private]

G4double BLCMDhelicaldipole::model [private]

G4double BLCMDhelicaldipole::psi0 [private]

G4double BLCMDhelicaldipole::phi0 [private]

G4double BLCMDhelicaldipole::Bsolenoid [private]

G4double BLCMDhelicaldipole::ez [private]

Referenced by construct(), and EndOfRunAction().

G4Polymarker BLCMDhelicaldipole::markers [private]

Referenced by BeginOfRunAction(), and EndOfRunAction().


The documentation for this class was generated from the following file:
g4beamline