BLCMDpillbox Class Reference

Inheritance diagram for BLCMDpillbox:

BLElement BLCommand

List of all members.


Detailed Description

BLCMDpillbox implements a pillbox.

Each placement of a BLCMDpillbox creates a PillboxField that is linked into BLGlobalField. Each PillboxField sets its timeOffset via the tune particle (re-stepping through its TimingVol if necessary).

If fieldMapFile is non-null, it is read as a BLFieldMap and determines the peak field; both B and E are multiplied by maxGradient of the cavity.

Public Member Functions

 BLCMDpillbox ()
 Default constructor. Defines the command, etc.
virtual ~BLCMDpillbox ()
 Destructor.
 BLCMDpillbox (const BLCMDpillbox &r)
 Copy constructor.
BLElementclone ()
 clone()
G4String commandName ()
 commandName() returns "pillbox".
int command (BLArgumentVector &argv, BLArgumentMap &namedArgs)
 command() implements the pillbox command.
void defineNamedArgs ()
 defineNamedArgs() defines the named arguments of the command.
void argChanged ()
 argChanged() does internal computations after some arg changed
void construct (G4RotationMatrix *relativeRotation, G4ThreeVector relativePosition, G4LogicalVolume *parent, G4String parentName, G4RotationMatrix *parentRotation, G4ThreeVector parentPosition)
 construct() - construct the pillbox. Creates a new PillboxField and adds it to BLGlobalField.
G4double getLength ()
 getLength() returns the overallLength of the pillbox
G4double getWidth ()
 getWidth() returns the outer radius of the pillbox
G4double getHeight ()
 getHeight() returns the outer radius of the pillbox
G4bool isOK ()
 isOK() returns true only if all placed pillboxes of this type have had their timing set, and if the tuning (if used) converged.
void generatePoints (int npoints, std::vector< G4ThreeVector > &v)
 generatePoints() from BLElement
G4bool isOutside (G4ThreeVector &local, G4double tolerance)
 isOutside() from BLElement

Private Attributes

G4double maxGradient
G4String color
G4double frequency
G4double innerLength
G4double innerRadius
G4double pipeThick
G4double wallThick
G4double irisRadius
G4double collarRadialThick
G4double collarThick
G4double win1Thick
G4double win1OuterRadius
G4double win2Thick
G4String winMat
G4double phaseAcc
G4double skinDepth
G4double timingTolerance
G4double maxStep
G4String cavityMaterial
G4double timeOffset
G4double timeIncrement
G4String fieldMapFile
G4double overallLength
G4double outerRadius
G4double omega
G4double rFactor
G4double E2Bfactor
BLFieldMapfieldMap
std::vector< PillboxField * > pillboxField

Friends

class PillboxField


Constructor & Destructor Documentation

BLCMDpillbox::BLCMDpillbox (  ) 

Default constructor. Defines the command, etc.

References BLCMDTYPE_ELEMENT, cavityMaterial, collarRadialThick, collarThick, color, E2Bfactor, fieldMap, fieldMapFile, frequency, innerLength, innerRadius, irisRadius, maxGradient, maxStep, omega, outerRadius, overallLength, phaseAcc, pipeThick, BLCommand::registerCommand(), rFactor, BLCommand::setDescription(), BLCommand::setSynopsis(), skinDepth, timeIncrement, timeOffset, timingTolerance, UNDETERMINED, wallThick, win1OuterRadius, win1Thick, win2Thick, and winMat.

Referenced by clone(), and command().

00186                            : BLElement(), pillboxField()
00187 {
00188         // register the commandName(), and its synopsis and description.
00189         registerCommand(BLCMDTYPE_ELEMENT);
00190         setSynopsis("Defines a pillbox RF cavity");
00191         setDescription("A Pillbox RF cavity is the basic RF element used to\n"
00192                 "construct a linac. The phaseAcc parameter sets the\n"
00193                 "phase of the tune particle at the center of the\n"
00194                 "cavity, and the timing offset of the cavity is determined\n"
00195                 "from that when the tune particle is tracked through\n"
00196                 "the cavity. Zero degrees is the rising zero-crossing of\n"
00197                 "the Ez field. If timeOffset is specified, it is used rather "
00198                 "than setting it from the Tune particle.\n\n"
00199                 "The Pipe, walls, and collars are always made of copper. "
00200                 "Pipe, wall, collar, win1, and win2 can be omitted by setting\n"
00201                 "their thickness to 0. Note that section 4.4 of the User's "
00202                 "Guide has a dimensioned drawing of a pillbox.");
00203 
00204         // provide initial values for fields
00205         maxGradient = UNDETERMINED; // will be converted to MV/m
00206         color = "1.0,0.0,0.0";
00207         frequency = UNDETERMINED;
00208         innerLength = UNDETERMINED;
00209         innerRadius = UNDETERMINED;
00210         pipeThick = 3.0*mm;
00211         wallThick = 3.0*mm;
00212         irisRadius = 11.0*cm;
00213         collarRadialThick = 5.0*mm;
00214         collarThick = 2.5*mm;
00215         win1Thick = 0.200*mm;
00216         win1OuterRadius = 5.0*cm;
00217         win2Thick = 0.500*mm;
00218         winMat = "Be";
00219         phaseAcc = 40.0*deg;
00220         skinDepth = 0.002*mm;
00221         timingTolerance = 0.001*ns;
00222         maxStep = -1.0;
00223         cavityMaterial = "Vacuum";
00224         timeOffset = UNDETERMINED;
00225         timeIncrement = 0.0;
00226         fieldMapFile = "";
00227         // non-arguments:
00228         overallLength = UNDETERMINED;
00229         outerRadius = UNDETERMINED;
00230         omega = UNDETERMINED;
00231         rFactor = UNDETERMINED;
00232         E2Bfactor = UNDETERMINED;
00233         fieldMap = 0;
00234 }

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

Destructor.

00133 { }

BLCMDpillbox::BLCMDpillbox ( const BLCMDpillbox r  ) 

Copy constructor.

References cavityMaterial, collarRadialThick, collarThick, color, BLTune::copyTunableArg(), E2Bfactor, fieldMap, fieldMapFile, frequency, innerLength, innerRadius, irisRadius, maxGradient, maxStep, omega, outerRadius, overallLength, phaseAcc, pipeThick, rFactor, skinDepth, timeIncrement, timeOffset, timingTolerance, wallThick, win1OuterRadius, win1Thick, win2Thick, and winMat.

00237                                                 : BLElement(r), pillboxField()
00238 {
00239         // copy fields one at a time (transfers default values from the
00240         // default object to this new object).
00241         BLTune::copyTunableArg(&maxGradient,&r.maxGradient);
00242         color = r.color;
00243         frequency = r.frequency;
00244         innerLength = r.innerLength;
00245         innerRadius = r.innerRadius;
00246         pipeThick = r.pipeThick;
00247         wallThick = r.wallThick;
00248         irisRadius = r.irisRadius;
00249         collarRadialThick = r.collarRadialThick;
00250         collarThick = r.collarThick;
00251         win1Thick = r.win1Thick;
00252         win1OuterRadius = r.win1OuterRadius;
00253         win2Thick = r.win2Thick;
00254         winMat = r.winMat;
00255         phaseAcc = r.phaseAcc;
00256         skinDepth = r.skinDepth;
00257         timingTolerance = r.timingTolerance;
00258         maxStep = r.maxStep;
00259         cavityMaterial = r.cavityMaterial;
00260         timeOffset = r.timeOffset;
00261         timeIncrement = r.timeIncrement;
00262         fieldMapFile = r.fieldMapFile;
00263 
00264         // non-arguments:
00265         overallLength = r.overallLength;
00266         outerRadius = r.outerRadius;
00267         omega = r.omega;
00268         rFactor = r.rFactor;
00269         E2Bfactor = r.E2Bfactor;
00270         fieldMap = r.fieldMap;
00271 }


Member Function Documentation

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

clone()

Implements BLElement.

References BLCMDpillbox().

00139 { return new BLCMDpillbox(*this); }

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

commandName() returns "pillbox".

Implements BLCommand.

00142 { return "pillbox"; }

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

command() implements the pillbox command.

Implements BLCommand.

References argChanged(), BLCMDpillbox(), cavityMaterial, E, E2Bfactor, fieldMap, fieldMapFile, frequency, g4bl_exit(), BLCommand::getMaterial(), BLCommand::handleNamedArgs(), innerRadius, maxGradient, BLCommand::print(), BLCommand::printError(), BLFieldMap::readFile(), and BLElement::setName().

00274 {
00275         // If the name of the cavity is "testUnits", perform the test and exit.
00276         if(argv[0] == "testUnits") {
00277                 // set the test case values, from Al Moretti's SuperFish example
00278                 maxGradient = 1.0 * megavolt/meter;
00279                 frequency = 0.956188/ns;
00280                 double Hcorrect = 1545.0; // A/m
00281                 double BcorrectTesla = 4.0*pi*1E-7 * Hcorrect; // MKS
00282                 argChanged();
00283                 // the first maximum of J1() has a value of 0.5819
00284                 double Bmax = maxGradient*E2Bfactor*0.5819;
00285                 printf("pillobx testUnits:\n           "         
00286                         "Freq(GHz)  Emax(MV/m)  Bmax(Tesla)  Radius(mm)\n");
00287                 printf("Values:  %10.4f %10.4f %12.4e %11.4f\n",
00288                         frequency/(1e9*hertz), maxGradient/(megavolt/meter),
00289                         Bmax/tesla, innerRadius/mm);
00290                 // these values already are in the correct units:
00291                 printf("Correct: %10.4f %10.4f %12.4e %8.1f\n",
00292                         0.956188,1.0,BcorrectTesla, 120.0);
00293                 printf("pillbox testUnits: exiting\n");
00294                 extern void g4bl_exit(int);
00295                 g4bl_exit(0);
00296         }
00297 
00298         if(argv.size() != 1) {
00299                 printError("pillbox: Invalid command, must have name");
00300                 return -1;
00301         }
00302 
00303         if(argv[0] == "default") {
00304                 return defaultPillbox.handleNamedArgs(namedArgs);
00305         }
00306 
00307         BLCMDpillbox *p = new BLCMDpillbox(defaultPillbox);
00308         p->setName(argv[0]);
00309         int retval = p->handleNamedArgs(namedArgs);
00310 
00311         if(p->fieldMapFile != "") {
00312                 p->fieldMap = new BLFieldMap();
00313                 p->fieldMap->readFile(p->fieldMapFile);
00314         }
00315 
00316         // check material exists
00317         if(p->cavityMaterial.size() > 0) getMaterial(p->cavityMaterial,false);
00318 
00319         p->print(argv[0]);
00320 
00321         return retval;
00322 }

void BLCMDpillbox::defineNamedArgs (  )  [virtual]

defineNamedArgs() defines the named arguments of the command.

Reimplemented from BLCommand.

References BLCommand::argDouble(), BLCommand::argString(), BLCommand::argTunable(), cavityMaterial, collarRadialThick, collarThick, color, fieldMapFile, frequency, innerLength, innerRadius, irisRadius, maxGradient, maxStep, phaseAcc, pipeThick, skinDepth, timeIncrement, timeOffset, timingTolerance, wallThick, win1OuterRadius, win1Thick, win2Thick, and winMat.

00325 {
00326         argTunable(maxGradient,"maxGradient","The peak gradient of the cavity (MV/m)",megavolt/meter);
00327         argString(color,"color","The color of the cavity");
00328         argDouble(frequency,"frequency","The frequency of the cavity (GHz)",1/ns,"",false);
00329         argDouble(innerLength,"innerLength","The inside length of the cavity (mm)",mm,"",false);
00330         argDouble(innerRadius,"innerRadius","The inside radius of the cavity (mm)",mm,"",false);
00331         argDouble(pipeThick,"pipeThick","The thickness of the pipe wall (mm)",mm,"",false);
00332         argDouble(wallThick,"wallThick","The thickness of the cavity walls (mm)",mm,"",false);
00333         argDouble(irisRadius,"irisRadius","The radius of the iris (mm)",mm,"",false);
00334         argDouble(collarRadialThick,"collarRadialThick","The radial thickness of the collar (mm)",mm,"",false);
00335         argDouble(collarThick,"collarThick","The thickness of the collar along z(mm)",mm,"",false);
00336         argDouble(win1Thick,"win1Thick","The thickness of the central portion of the\n"
00337                 "windows; zero for no window (mm)",mm,"",false);
00338         argDouble(win1OuterRadius,"win1OuterRadius","The radius of the central portion of\n"
00339                 "the windows (mm)",mm,"",false);
00340         argDouble(win2Thick,"win2Thick","The thickness of the outer portion of the\n"
00341                 "windows; zero for no window (mm)",mm,"",false);
00342         argString(winMat,"winMat","The material of the windows");
00343         argDouble(phaseAcc,"phaseAcc","The reference phase of the cavity (degrees)",deg);
00344         argDouble(skinDepth,"skinDepth","The skin depth (mm)",mm,"",false);
00345         argDouble(timingTolerance,"timingTolerance","Tolerance for timing tuning (ns)",ns);
00346         argDouble(maxStep,"maxStep","The maximum stepsize in the element (mm).",mm);
00347         argString(cavityMaterial,"cavityMaterial","Material of cavity volume (Vacuum).");
00348         argDouble(timeOffset,"timeOffset","Time offset for cavity (default: tuned by tune particle) (ns).",ns);
00349         argDouble(timeIncrement,"timeIncrement","Increment to timeOffset, applied AFTER tuning. (ns).",ns);
00350         argString(fieldMapFile,"fieldMapFile","Filename for BLFieldMap (pillbox if null).",false);
00351 }

void BLCMDpillbox::argChanged (  )  [virtual]

argChanged() does internal computations after some arg changed

Reimplemented from BLCommand.

References collarThick, E2Bfactor, frequency, BLParam::getDouble(), innerLength, innerRadius, maxStep, omega, outerRadius, overallLength, Param, pipeThick, BLCommand::printError(), rFactor, and UNDETERMINED.

Referenced by command().

00354 {
00355         // compute non-argument variables
00356         omega = 2.0 * pi * frequency;
00357         if(innerRadius == UNDETERMINED) {
00358                 if(frequency > 0.0) {
00359                         // the first zero of J0(r) is at r=2.405
00360                         innerRadius =  2.405 / (omega/c_light);
00361                 } else {
00362                         printError("pillbox: innerRadius is undefined");
00363                 }
00364         }
00365         overallLength = innerLength + 2.0 * collarThick;
00366         outerRadius = innerRadius + pipeThick;
00367         rFactor = omega/c_light;
00368         // This is the incorrect value from the BeamTools:
00369         //     E2Bfactor = 1.0e-9 / c_light;   // MeV/m, MKS units
00370         // C.f. Eq. 8.77 in J.D.Jackson, _Classical_Electrodynamics_,
00371         //       first edition (1962). His value is 1/c in vacuum.
00372         E2Bfactor = 1.0/c_light;   // Geant4 units
00373         if(maxStep < 0.0) maxStep = Param.getDouble("maxStep");
00374 }

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

construct() - construct the pillbox. Creates a new PillboxField and adds it to BLGlobalField.

Implements BLElement.

References BLGlobalField::addElementField(), cavityMaterial, collarRadialThick, collarThick, BLCommand::getMaterial(), BLElement::getName(), BLManager::getObject(), BLGlobalField::getObject(), BLCoordinateTransform::getPosition(), PillboxField::global2local, innerLength, innerRadius, irisRadius, BLCoordinateTransform::isRotated(), maxStep, PillboxField::pillbox, pillboxField, PillboxField, pipeThick, BLManager::registerTuneParticleStep(), PillboxField::timingPV, wallThick, win1OuterRadius, win1Thick, win2Thick, and winMat.

00382 {
00383         G4String thisname = parentName+getName();
00384 
00385         G4Material *cu = getMaterial("Cu");
00386         G4Material *wm = getMaterial(winMat);
00387         G4Material *cavityMat = getMaterial(cavityMaterial);
00388         G4UserLimits *userLimits = new G4UserLimits(maxStep);
00389 
00390         // find enclosing Tubs
00391         G4double halflen = innerLength/2+wallThick;
00392         halflen = std::max(halflen,innerLength/2+collarThick);
00393         halflen = std::max(halflen,innerLength/2+win1Thick);
00394         halflen = std::max(halflen,innerLength/2+win2Thick);
00395         G4Tubs *tubs = new G4Tubs(getName(),0.0,innerRadius+pipeThick,
00396                                         halflen,0.0,2.0*pi);
00397         G4LogicalVolume *pillbox = new G4LogicalVolume(tubs,cavityMat,
00398                                         parentName+getName(),0,0,userLimits);
00399 
00400         G4ThreeVector loc;
00401 
00402         // The outer pipe
00403         G4LogicalVolume *pipe=0;
00404         if(pipeThick > 0.001*mm) {
00405                 tubs = new G4Tubs(getName()+"pipe",innerRadius,
00406                                         innerRadius+pipeThick,
00407                                         innerLength/2+wallThick,0.0,2.0*pi);
00408                 pipe = new G4LogicalVolume(tubs,cu,
00409                                         getName()+"pipeLog",0,0,userLimits);
00410                 new G4PVPlacement(0,loc,pipe,parentName+getName()+"Pipe",
00411                                         pillbox,false,0);
00412         }
00413 
00414         // the wall (us and ds) - radially: collarOuterRadius to innerRadius
00415         G4LogicalVolume *wall=0;
00416         if(wallThick > 0.001*mm) {
00417                 tubs = new G4Tubs(getName()+"wall",irisRadius+collarRadialThick,
00418                                         innerRadius,wallThick/2,0.0,2.0*pi);
00419                 wall = new G4LogicalVolume(tubs,cu,
00420                                         getName()+"wallLog",0,0,userLimits);
00421                 G4double wallCenterZ = innerLength/2.0 + wallThick/2.0;
00422                 loc.setZ(-wallCenterZ);
00423                 new G4PVPlacement(0,loc,wall,parentName+getName()+"UsWall",
00424                                         pillbox,false,0);
00425                 loc.setZ(+wallCenterZ);
00426                 new G4PVPlacement(0,loc,wall,parentName+getName()+"DsWall",
00427                                         pillbox,false,0);
00428         }
00429         
00430         // the collar (us and ds) - includes wall thickness (i.e. the inner
00431         // edge of the collar is the same as the inner edge of the wall)
00432         G4LogicalVolume *collar=0;
00433         if(collarThick > 0.001*mm) {
00434                 tubs = new G4Tubs(getName()+"collar",irisRadius,
00435                                         irisRadius+collarRadialThick,
00436                                         collarThick/2,0.0,2.0*pi);
00437                 collar = new G4LogicalVolume(tubs,cu,
00438                                         getName()+"collarLog",0,0,userLimits);
00439                 G4double collarCenterZ = innerLength/2.0 + collarThick/2.0;
00440                 loc.setZ(-collarCenterZ);
00441                 new G4PVPlacement(0,loc,collar,parentName+getName()+"UsCollar",
00442                                         pillbox,false,0);
00443                 loc.setZ(+collarCenterZ);
00444                 new G4PVPlacement(0,loc,collar,parentName+getName()+"DsCollar",
00445                                         pillbox,false,0);
00446         }
00447         
00448         // win1 (the inner window)
00449         G4LogicalVolume *win1=0;
00450         if(win1Thick > 0.001*mm) {
00451                 tubs = new G4Tubs(getName()+"win1",0.0,win1OuterRadius,
00452                                         win1Thick/2,0.0,2.0*pi);
00453                 win1 = new G4LogicalVolume(tubs,wm,
00454                                         getName()+"win1Log",0,0,userLimits);
00455                 G4double win1CenterZ = innerLength/2.0 + win1Thick/2.0;
00456                 loc.setZ(-win1CenterZ);
00457                 new G4PVPlacement(0,loc,win1,parentName+getName()+"UsWin1",
00458                                         pillbox,false,0);
00459                 loc.setZ(+win1CenterZ);
00460                 new G4PVPlacement(0,loc,win1,parentName+getName()+"DsWin1",
00461                                         pillbox,false,0);
00462         }
00463 
00464         // win2 (the outer window)
00465         G4LogicalVolume *win2=0;
00466         if(win2Thick > 0.001*mm) {
00467                 tubs = new G4Tubs(getName()+"win2",win1OuterRadius,irisRadius,
00468                                         win2Thick/2,0.0,2.0*pi);
00469                 win2 = new G4LogicalVolume(tubs,wm,
00470                                         getName()+"win2Log",0,0,userLimits);
00471                 G4double win2CenterZ = innerLength/2.0 + win2Thick/2.0;
00472                 loc.setZ(-win2CenterZ);
00473                 new G4PVPlacement(0,loc,win2,parentName+getName()+"UsWin2",
00474                                         pillbox,false,0);
00475                 loc.setZ(+win2CenterZ);
00476                 new G4PVPlacement(0,loc,win2,parentName+getName()+"DsWin2",
00477                                         pillbox,false,0);
00478         }
00479 
00480         // half of the interior
00481         tubs = new G4Tubs(getName()+"half",0.0,innerRadius,
00482                                         innerLength/4,0.0,2.0*pi);
00483         G4LogicalVolume *half = new G4LogicalVolume(tubs,cavityMat,getName()+"halfLog",
00484                                         0,0,userLimits);
00485         G4double halfCenterZ = innerLength/4.0;
00486         loc.setZ(+halfCenterZ);
00487         new G4PVPlacement(0,loc,half,parentName+getName()+"DsHalf",
00488                                                         pillbox,false,0);
00489         loc.setZ(-halfCenterZ);
00490         G4VPhysicalVolume *timingPV = new G4PVPlacement(0,loc,half,
00491                                 parentName+getName()+"TimingVol",pillbox,false,0);
00492 
00493 #ifdef G4VIS_USE
00494         G4VisAttributes *c1 = new G4VisAttributes(true,G4Color(0.5,0.0,0.0));
00495         G4VisAttributes *c2 = new G4VisAttributes(true,G4Color(0.7,0.0,0.0));
00496         G4VisAttributes *c3 = new G4VisAttributes(true,G4Color(1.0,0.0,0.0));
00497         const G4VisAttributes *invisible = &G4VisAttributes::GetInvisible();
00498         if(pipe) pipe->SetVisAttributes(c1);
00499         if(wall) wall->SetVisAttributes(c1);
00500         if(collar) collar->SetVisAttributes(c1);
00501         if(win1) win1->SetVisAttributes(c3);
00502         if(win2) win2->SetVisAttributes(c2);
00503         half->SetVisAttributes(invisible);
00504         pillbox->SetVisAttributes(invisible);
00505 #endif
00506 
00507         // geant4 rotation convention is backwards from g4beamline
00508         G4RotationMatrix *g4rot = 0;
00509         if(relativeRotation)
00510                 g4rot = new G4RotationMatrix(relativeRotation->inverse());
00511 
00512         // place the pillbox into the parent
00513         new G4PVPlacement(g4rot,relativePosition,pillbox,
00514                                 parentName+getName(),parent,false,0);
00515 
00516         // get globalRotation and globalPosition
00517         G4RotationMatrix *globalRotation = 0;
00518         if(relativeRotation && parentRotation) {
00519                 globalRotation = 
00520                     new G4RotationMatrix(*parentRotation * *relativeRotation);
00521         } else if(relativeRotation) {
00522                 globalRotation = relativeRotation;
00523         } else if(parentRotation) {
00524                 globalRotation = parentRotation;
00525         }
00526         G4ThreeVector globalPosition(relativePosition + parentPosition);
00527         if(parentRotation)
00528                 globalPosition = *parentRotation * relativePosition +
00529                                 parentPosition;
00530 
00531         BLCoordinateTransform global2local(globalRotation,globalPosition);
00532 
00533         // add this pillbox field to the GlobalField, and to BLManager
00534         G4double dz = innerLength/2.0 + wallThick;
00535         // if rotated, make an overestimate of the field occupancy along z
00536         if(global2local.isRotated())
00537                 dz += innerRadius;
00538         PillboxField *pf = 
00539                 new PillboxField(thisname,global2local,this,timingPV);
00540         BLGlobalField::getObject()->addElementField(pf);
00541         pillboxField.push_back(pf);
00542         BLManager::getObject()->registerTuneParticleStep(timingPV,pf);
00543 
00544         printf("BLCMDpillbox::construct %s parent=%s relZ=%.1f globZ=%.1f\n"
00545                         "\tzmin=%.1f zmax=%.1f\n",
00546                 thisname.c_str(),parentName.c_str(),relativePosition[2],
00547                 global2local.getPosition()[2], global2local.getPosition()[2]-dz,
00548                 global2local.getPosition()[2]+dz);
00549 }

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

getLength() returns the overallLength of the pillbox

Implements BLElement.

References overallLength.

00163 { return overallLength; }

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

getWidth() returns the outer radius of the pillbox

Implements BLElement.

References outerRadius.

00166 { return outerRadius*2.0; }

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

getHeight() returns the outer radius of the pillbox

Implements BLElement.

References outerRadius.

00169 { return outerRadius*2.0; }

G4bool BLCMDpillbox::isOK (  )  [virtual]

isOK() returns true only if all placed pillboxes of this type have had their timing set, and if the tuning (if used) converged.

Implements BLElement.

References isUndetermined(), and pillboxField.

00552 {
00553         G4bool retval = true;
00554 
00555         // verify all placed pillboxField-s have set their timeOffset 
00556         std::vector<PillboxField *>::iterator i;
00557         for(i=pillboxField.begin(); i!=pillboxField.end(); ++i) {
00558                 if(isUndetermined((*i)->timeOffset) || (*i)->settingTiming) {
00559                         printf("*** BLCMDpillbox %s has not determined timing\n",
00560                                 (*i)->getName().c_str());
00561                         retval = false;
00562                 }
00563         }
00564 
00565         return retval;
00566 }

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

generatePoints() from BLElement

Implements BLElement.

References BLElement::generateTubs(), outerRadius, and overallLength.

00759 {
00760         generateTubs(npoints, 0.0, outerRadius, 0.0, 360.0*deg,
00761                         overallLength, v);
00762 }

G4bool BLCMDpillbox::isOutside ( G4ThreeVector &  local,
G4double  tolerance 
) [virtual]

isOutside() from BLElement

Implements BLElement.

References outerRadius, and overallLength.

00765 {
00766         G4double r = sqrt(local[0]*local[0]+local[1]*local[1]);
00767         return r > outerRadius-tolerance ||
00768                 fabs(local[2]) > overallLength/2.0-tolerance;
00769 }


Friends And Related Function Documentation

friend class PillboxField [friend]

Referenced by construct().


Member Data Documentation

G4double BLCMDpillbox::maxGradient [private]

G4String BLCMDpillbox::color [private]

Referenced by BLCMDpillbox(), and defineNamedArgs().

G4double BLCMDpillbox::frequency [private]

G4double BLCMDpillbox::innerLength [private]

G4double BLCMDpillbox::innerRadius [private]

G4double BLCMDpillbox::pipeThick [private]

G4double BLCMDpillbox::wallThick [private]

G4double BLCMDpillbox::irisRadius [private]

G4double BLCMDpillbox::collarRadialThick [private]

G4double BLCMDpillbox::collarThick [private]

G4double BLCMDpillbox::win1Thick [private]

G4double BLCMDpillbox::win1OuterRadius [private]

G4double BLCMDpillbox::win2Thick [private]

G4String BLCMDpillbox::winMat [private]

G4double BLCMDpillbox::phaseAcc [private]

G4double BLCMDpillbox::skinDepth [private]

G4double BLCMDpillbox::timingTolerance [private]

G4double BLCMDpillbox::maxStep [private]

G4String BLCMDpillbox::cavityMaterial [private]

G4double BLCMDpillbox::timeOffset [private]

G4double BLCMDpillbox::timeIncrement [private]

G4String BLCMDpillbox::fieldMapFile [private]

G4double BLCMDpillbox::overallLength [private]

G4double BLCMDpillbox::outerRadius [private]

G4double BLCMDpillbox::omega [private]

G4double BLCMDpillbox::rFactor [private]

G4double BLCMDpillbox::E2Bfactor [private]

std::vector<PillboxField *> BLCMDpillbox::pillboxField [private]

Referenced by construct(), and isOK().


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