A cornerarc is an approximation to an arc in the centerline, implemented by 3 corners such that the overall bend angle and centerline path length are the same as the arc would be.
The error (largest distance between the arc and the actual centerline) increases with angle, but is less than 10% of radius for 90 degrees.
Public Member Functions | |
BLCMDcornerarc () | |
Default constructor. | |
~BLCMDcornerarc () | |
Destructor. | |
BLCMDcornerarc (const BLCMDcornerarc &r) | |
Copy constructor. | |
G4String | commandName () |
commandName() returns "cornerarc". | |
int | command (BLArgumentVector &argv, BLArgumentMap &namedArgs) |
command() implements the cornerarc command. | |
void | defineNamedArgs () |
defineNamedArgs() defines the named arguments for the command. | |
Private Attributes | |
G4double | z |
G4double | centerRadius |
G4double | angle |
G4String | rotation |
G4double | radiusCut |
G4double | finalZ |
BLCMDcornerarc::BLCMDcornerarc | ( | ) |
Default constructor.
References angle, BLCMDTYPE_LAYOUT, centerRadius, finalZ, radiusCut, BLCommand::registerCommand(), rotation, BLCommand::setDescription(), BLCommand::setSynopsis(), and z.
Referenced by command().
00076 : BLCommand(), rotation() 00077 { 00078 registerCommand(BLCMDTYPE_LAYOUT); 00079 setSynopsis("Implement a cornerarc in the centerline."); 00080 setDescription("The centerline is bent by a rotation.\n" 00081 "Three corners are used to approximate an arc; the total\n" 00082 "angle and path-length are equal to those for the arc.\n" 00083 "Should be used immediately after an idealsectorbend or\n" 00084 "a genericbend.\n" 00085 "The z value is for the front face of the arc.\n" 00086 "\nNOTE: This command is self-placing, do not use the place\n" 00087 "command; it also affects all following placements, and\n" 00088 "it cannot be issued inside a group.\n\n" 00089 "This command is well matched to a sector bend, but can also " 00090 "be used with a normal bending magnet -- normally the magnet " 00091 "is placed before the cornerarc and is rotated by half the " 00092 "bend angle.\n\n" 00093 "The only useful rotations are those around the centerline Z.\n" 00094 "\nNOTE: all placements before this command must have z\n" 00095 "values before the corner, and all placements after this\n" 00096 "command must have z values after the corner. Note also\n" 00097 "that the angle is limited to 90 degrees.\n\n" 00098 "Note that the radiusCut is important to reduce or eliminate " 00099 "ambiguities in the global to centerline coordinate transform. " 00100 "It can also be used to 'shield' the beamline to prevent " 00101 "particles from taking unusual paths around the outside " 00102 "of beamline elements."); 00103 // initial field values 00104 z = 0.0; 00105 centerRadius=0.0; 00106 angle = 0.0; 00107 rotation = ""; 00108 radiusCut = 0.0; 00109 finalZ = 0.0; 00110 }
BLCMDcornerarc::BLCMDcornerarc | ( | const BLCMDcornerarc & | r | ) |
G4String BLCMDcornerarc::commandName | ( | ) | [inline, virtual] |
int BLCMDcornerarc::command | ( | BLArgumentVector & | argv, | |
BLArgumentMap & | namedArgs | |||
) | [virtual] |
command() implements the cornerarc command.
Implements BLCommand.
References angle, BLCMDcornerarc(), centerRadius, BLCoordinates::corner(), finalZ, BLGroup::getCurrent(), BLGroup::getWorld(), BLCommand::handleNamedArgs(), BLCommand::print(), BLCommand::printError(), radiusCut, rotation, BLCommand::stringToRotationMatrix(), and z.
00126 { 00127 if(argv.size() >= 1 && argv[0] == "default") { 00128 return defaultCornerArc.handleNamedArgs(namedArgs); 00129 } 00130 00131 if(BLGroup::getWorld() != BLGroup::getCurrent()) { 00132 printError("cornerarc: cannot be issued inside a group."); 00133 return -1; 00134 } 00135 00136 BLCMDcornerarc *t = new BLCMDcornerarc(defaultCornerArc); 00137 t->z = -DBL_MAX; 00138 int retval = t->handleNamedArgs(namedArgs); 00139 00140 if(t->z == -DBL_MAX) { 00141 printError("cornerarc: z must be specified"); 00142 return -1; 00143 } 00144 00145 G4double angleAbs = fabs(t->angle); 00146 00147 if(angleAbs < 0.001 || angleAbs > 3.1416/2.0) { 00148 printError("cornerarc: invalid angle (0<|angle|<=90 degrees)"); 00149 return -1; 00150 } 00151 if(t->centerRadius < 0.001*mm) { 00152 printError("cornerarc: invalid centerRadius"); 00153 return -1; 00154 } 00155 00156 // Compute the angles for the 3 corners that approximate the arc 00157 G4double beta = acos(sin(angleAbs/2.0)/(angleAbs/2.0)); 00158 G4double phi1 = angleAbs/2.0 - beta; 00159 G4double phi2 = beta + beta; 00160 G4double deltaz = t->centerRadius * angleAbs / 2.0; 00161 00162 // get the rotationMatrix-s 00163 char temp[64]; 00164 int i = sprintf(temp,"Y%.6f",(t->angle>0.0 ? phi1 : -phi1)/deg); 00165 if(t->rotation.size() > 0) 00166 sprintf(temp+i,",%s",t->rotation.c_str()); 00167 G4RotationMatrix matrix1 = *stringToRotationMatrix(temp); 00168 i = sprintf(temp,"Y%.6f",(t->angle>0.0 ? phi2 : -phi2)/deg); 00169 G4RotationMatrix matrix2 = *stringToRotationMatrix(temp); 00170 i = sprintf(temp,"Y%.6f",(t->angle>0.0 ? phi1 : -phi1)/deg); 00171 G4RotationMatrix matrix3 = *stringToRotationMatrix(temp); 00172 00173 BLCoordinates::corner(t->z,&matrix1,t->radiusCut); 00174 BLCoordinates::corner(t->z+deltaz,&matrix2,t->radiusCut); 00175 BLCoordinates::corner(t->z+deltaz+deltaz,&matrix3,t->radiusCut); 00176 00177 t->finalZ = t->z+deltaz+deltaz; 00178 00179 t->print(""); 00180 00181 return retval; 00182 }
void BLCMDcornerarc::defineNamedArgs | ( | ) | [virtual] |
defineNamedArgs() defines the named arguments for the command.
Reimplemented from BLCommand.
References angle, BLCommand::argDouble(), BLCommand::argString(), centerRadius, radiusCut, rotation, and z.
00185 { 00186 argDouble(z,"z","The centerline Z of the cornerarc (mm).",mm); 00187 argDouble(centerRadius,"centerRadius","The radius of the centerline arc (mm).",mm); 00188 argDouble(angle,"angle","The angle of bend, >0=left, <0=right (degrees).",deg); 00189 argString(rotation,"rotation","The rotation of the cornerarc (see above)."); 00190 argDouble(radiusCut,"radiusCut","The radius cut for this following segment (mm).",mm); 00191 }
G4double BLCMDcornerarc::z [private] |
Referenced by BLCMDcornerarc(), command(), and defineNamedArgs().
G4double BLCMDcornerarc::centerRadius [private] |
Referenced by BLCMDcornerarc(), command(), and defineNamedArgs().
G4double BLCMDcornerarc::angle [private] |
Referenced by BLCMDcornerarc(), command(), and defineNamedArgs().
G4String BLCMDcornerarc::rotation [private] |
Referenced by BLCMDcornerarc(), command(), and defineNamedArgs().
G4double BLCMDcornerarc::radiusCut [private] |
Referenced by BLCMDcornerarc(), command(), and defineNamedArgs().
G4double BLCMDcornerarc::finalZ [private] |
Referenced by BLCMDcornerarc(), and command().