DECLARE SUB dpoldatbl (datfile$, bl!, sbl!, th!, sth!)
DECLARE SUB dpolgetblav (bl!, sbl!, th!, sth!)
DECLARE SUB dpolpltbl (pltfile$, bl!, sbl!, th!, sth!)
DECLARE SUB bllogcalc (logfile$)
DECLARE SUB blngetcoilradius (rcoil!)
DECLARE SUB imaggetiav (imag!, simag!)
DECLARE SUB imagramp (imagtest!)
DECLARE SUB muxsetchan (chan$)
DECLARE SUB vtautogain ()
DECLARE SUB vtgetvthar2 (rotdir$, vtmag!, vtphase!)
DECLARE SUB qpolcalcxy (nh%, bl!(), sbl!(), th!(), sth!(), x!, sx!, y!, sy!)
DECLARE SUB qpolcheckpar (ok$)
DECLARE SUB qpoldatgl (datfile$, gl!, sgl!, th!, sth!)
DECLARE SUB qpoldatglvsi (datfile$, imag!, simag!, gl!, sgl!, tf!, stf!)
DECLARE SUB qpoldatxy (logfile$, th!(), x!, sx!, y!, sy!)
DECLARE SUB qpolgetgl (gl!, th!)
DECLARE SUB qpolgetglav (gl!, sgl!, th!, sth!)
DECLARE SUB qpollogcalc (logfile$)
DECLARE SUB qpolloggl (logfile$, coilconst!, vt2!, vtphase!, gl!, th!)
DECLARE SUB qpollogglav (logfile$, gl!, sgl!, th!, sth!)
DECLARE SUB qpollogstrmeas (logfile$)
DECLARE SUB qpolpltgl (pltfile$, gl!, sgl!, th!, sth!)
DECLARE SUB qpolpltglvsi (pltfile$, imag!, gl!, sgl!, tf!, stf!)
DECLARE SUB blcalcratio (nh%, bln!(), sbln!(), nmain%, rn!(), srn!())
DECLARE SUB bldatblhar (datfile$, nmain%, nh%, bln!(), sbln!(), thn!(), sthn!(), r!(), sr!())
DECLARE SUB bldatblvsihar (datfile$, imag!, simag!, nmain%, nh%, bln!(), sbln!(), thn!(), sthn!(), r!(), sr!())
DECLARE SUB blgetblnav (coilconfig$, nh%, bln!(), sbln!(), thn!(), sthn!())
DECLARE SUB blpltblhar (pltfile$, nh%, rn!(), srn!())
DECLARE SUB blpltblvsihar (pltfile$, legend$, nh%, rn!(), srn!())
DECLARE SUB pmtdattslice (logfile$, name1$, t1!, name2$, t2!, name3$, t3!)
DECLARE SUB pmtgettslice (pmt1nameP$, t1!, pmt2nameP$, t2!, pmt3nameP$, t3!)
DECLARE SUB filenewpage (filename$)

'****************************************************************************
'Module QPOL.BAS
'
'These subroutines perform calculations related to quadrupole magnet
'measurements.
'
'Zachary Wolf
'1/2/95
'****************************************************************************

'Put the required parameters for this module in a common block
COMMON SHARED /qpol/ logfileP$, strdatfileP$, strpltfileP$, hardatfileP$, harpltfileP$, magnameP$, runP$, strqpolcoilconstP!, blnmeasaveP%, nimagtestP%, imagtestP!(), imagharflagP!(), nhardisplayP%

SUB qpolcalcxy (nh%, bl!(), sbl!(), th!(), sth!(), x!, sx!, y!, sy!)

'****************************************************************************
'This subroutine calculates the (x, y) position of the quadrupole magnetic
'center relative to the coil axis.  X is horizontal to the right (viewed
'from the encoder end) and Y is up.
'
'Input:
' nh%, the number of voltage harmonics recorded
' bl!(1 to nh%), the integrated strength of each harmonic at the coil radius (T)
' sbl!(1 to nh%), rms variation of bl
' th!(1 to nh%), the angle of the first south pole of each harmonic wrt horizontal (deg)
' sth!(1 to nh%), rms variation of th
'
'Output:
' x!, x coordinate of the magnetic center wrt the coil axis (m)
' sx!, calculated error for x
' y!, y coordinate of the magnetic center wrt the coil axis (m)
' sy!, calculated error for y
'
'Zachary Wolf
'1/2/94
'****************************************************************************

'Simplify the notation
pi! = 3.1415926#

'Coil parameters
CALL blngetcoilradius(rcoil!)

'Compute gl from bl(2)
gl! = bl!(2) / rcoil!
sgl! = sbl!(2) / rcoil!

'Compute the magnetic center relative to the coil
'For a normal quad (th(2) = 45 deg):
x! = -(1! / gl!) * bl!(1) * SIN(th!(1) * pi! / 180!)
y! = -(1! / gl!) * bl!(1) * COS(th!(1) * pi! / 180!)
'With the leads reversed, the signs change
IF th!(2) < 0 THEN
  x! = -x!
  y! = -y!
END IF

'Compute the error on the magnetic center
sx! = ABS(x!) * SQR((sgl! / gl!) ^ 2 + (sbl!(1) / bl!(1)) ^ 2 + (COS(th!(1) * pi! / 180!) / (SIN(th!(1) * pi! / 180!) + .000001)) ^ 2 * (sth!(1) * pi! / 180!) ^ 2)
sy! = ABS(y!) * SQR((sgl! / gl!) ^ 2 + (sbl!(1) / bl!(1)) ^ 2 + (SIN(th!(1) * pi! / 180!) / (COS(th!(1) * pi! / 180!) + .000001)) ^ 2 * (sth!(1) * pi! / 180!) ^ 2)

END SUB

SUB qpolcheckpar (ok$)

'****************************************************************************
'This subroutine checks that all required parameters have been assigned
'values.
'
'Output:
' ok$, "y" if parameter assignments are ok, "n" otherwise
'
'Zachary Wolf
'9/18/95
'****************************************************************************

'Default value
ok$ = "n"

'Files
IF logfileP$ = "" OR strdatfileP$ = "" OR strpltfileP$ = "" OR hardatfileP$ = "" OR harpltfileP$ = "" THEN
 PRINT "QPOL: required file not defined"
 EXIT SUB
END IF

'Runinfo
IF magnameP$ = "" OR runP$ = "" THEN
 PRINT "QPOL: required run parameter not defined"
 EXIT SUB
END IF

'Coil
IF strqpolcoilconstP! <= 0 THEN
 PRINT "QPOL: required coil parameter not defined"
 EXIT SUB
END IF

'Measurement parameters
IF blnmeasaveP% <= 0 OR nhardisplayP% <= 0 OR nimagtestP% <= 0 THEN
 PRINT "QPOL: required measurement parameter not defined"
 EXIT SUB
END IF

'If we made it this far, all parameters have values
ok$ = "y"

END SUB

SUB qpoldatgl (datfile$, gl!, sgl!, th!, sth!)

'****************************************************************************
'This subroutine writes a gl measurement to the data file.
'
'Input:
' datfile$, the name of the data file
' gl!, the integrated gradient at the measured current
' sgl!, the standard deviation of gl
' th!, south pole angle (deg)
' sth!, rms variation of TH
'
'Zachary Wolf
'11/19/95, 8/19/96, 10/23/96
'****************************************************************************

'Keep track of Ncall so a header can be written on the first call
STATIC ncall%
ncall% = ncall% + 1

'Open the log file
filenum% = FREEFILE
OPEN datfile$ FOR APPEND AS filenum%

'On the first call, write a header for the hysteresis loop
IF ncall% = 1 THEN
' PRINT #filenum%,
' PRINT #filenum%,
' PRINT #filenum%, " Coil Constant = "; strqpolcoilconstP!; " 1/m^2"
 PRINT #filenum%,
 PRINT #filenum%,
 PRINT #filenum%, "  INTEGRATED GRADIENT AND SOUTH POLE ANGLE  "
 PRINT #filenum%,
 PRINT #filenum%, "     GL          sigGL       TH       sTH   "
 PRINT #filenum%, "     (T)          (T)       (deg)    (deg)  "
 PRINT #filenum%, "------------+------------  --------+--------"
END IF

'Write the result
PRINT #filenum%, USING "#.######^^^^"; gl!;
PRINT #filenum%, USING " #.######^^^^"; sgl!;
PRINT #filenum%, USING "  ####.###"; th!;
PRINT #filenum%, USING " ####.###"; sth!

'Close the log file
CLOSE filenum%

END SUB

SUB qpoldatglvsi (datfile$, imag!, simag!, gl!, sgl!, tf!, stf!)

'****************************************************************************
'This subroutine writes a gl measurement at Imag to the data file.
'
'Input:
' datfile$, the name of the data file
' imag!, the measured current
' simag!, the standard deviation of imag
' gl!, the integrated gradient at the measured current
' sgl!, the standard deviation of gl
' tf!, GL/I
' stf!, the standard deviation of TF
'
'Zachary Wolf
'11/19/95
'****************************************************************************

'Keep track of Ncall so a header can be written on the first call
STATIC ncall%
ncall% = ncall% + 1

'Open the log file
filenum% = FREEFILE
OPEN datfile$ FOR APPEND AS filenum%

'On the first call, write a header for the hysteresis loop
IF ncall% = 1 THEN
 PRINT #filenum%,
 PRINT #filenum%,
 PRINT #filenum%,
 PRINT #filenum%, "             INTEGRATED GRADIENT VS CURRENT"
 PRINT #filenum%,
 PRINT #filenum%, "  Imag    sigImag     GL     sigGL      GL/I   sigGL/I "
 PRINT #filenum%, "  (A)       (A)       (T)     (T)      (T/kA)   (T/kA) "
 PRINT #filenum%, "--------+--------  --------+--------  --------+--------"
END IF

'Write the result
PRINT #filenum%, USING "####.###"; imag!;
PRINT #filenum%, USING " ####.###"; simag!;
PRINT #filenum%, USING "  ##.#####"; gl!;
PRINT #filenum%, USING " ##.#####"; sgl!;
PRINT #filenum%, USING "  ###.####"; tf! * 1000!;
PRINT #filenum%, USING " ###.####"; stf! * 1000!

'Close the log file
CLOSE filenum%

END SUB

SUB qpoldatxy (logfile$, th!(), x!, sx!, y!, sy!)

'****************************************************************************
'This subroutine writes the (x, y) location of the magnetic center
'relative to the coil to the log file.
'
'Input:
' logfile$, file to print results
' th!(1 to nh%), gives the orientation of the quadrupole field
' x!, x position of the magnetic center wrt the coil axis
' sx!, error on x
' y!, y position of the magnetic center wrt the coil axis
' sy!, error on y
'
'Zachary Wolf
'1/4/95
'****************************************************************************

'Open the log file
filenum% = FREEFILE
OPEN logfile$ FOR APPEND AS filenum%

'Print the magnetic center to the log file
PRINT #filenum%,
PRINT #filenum%,
PRINT #filenum%, "The magnetic center is at"
PRINT #filenum%, "X = ";
PRINT #filenum%, USING "######.#"; x! * 1000000!;
PRINT #filenum%, " +- ";
PRINT #filenum%, USING "######.#"; sx! * 1000000!;
PRINT #filenum%, " microns"
PRINT #filenum%, "Y = ";
PRINT #filenum%, USING "######.#"; y! * 1000000!;
PRINT #filenum%, " +- ";
PRINT #filenum%, USING "######.#"; sy! * 1000000!;
PRINT #filenum%, " microns"
'PRINT #filenum%,
'PRINT #filenum%, "The magnetic center is measured relative to the coil axis"
'PRINT #filenum%, "effectively at the longitudinal center of the magnet."
'PRINT #filenum%, "Pitch and yaw are not seen by the coil."

'Print the following diagram
PRINT #filenum%,
PRINT #filenum%, "View from the non-lead end of the magnet:"
PRINT #filenum%,
IF th!(2) > 0 THEN
  PRINT #filenum%, "           y          "
  PRINT #filenum%, "           ^          "
  PRINT #filenum%, "           |          "
  PRINT #filenum%, "        N     S       "
  PRINT #filenum%,
  PRINT #filenum%, "                --> x "
  PRINT #filenum%,
  PRINT #filenum%, "        S     N       "
ELSE
  PRINT #filenum%, "           y          "
  PRINT #filenum%, "           ^          "
  PRINT #filenum%, "           |          "
  PRINT #filenum%, "        S     N       "
  PRINT #filenum%,
  PRINT #filenum%, "                --> x "
  PRINT #filenum%,
  PRINT #filenum%, "        N     S       "
END IF

'Close the log file
CLOSE filenum%

END SUB

SUB qpolgetgl (gl!, th!)

'****************************************************************************
'This subroutine calculates the integrated quadrupole strength from the
'measured coil voltage at the second harmonic.  It also calculates the
'angle of the first south pole.
'
'Output:
' gl!, the integrated gradient (T)
' th!, the south pole angle (deg)
'
'Zachary Wolf
'11/19/95
'****************************************************************************

'Make sure all parameters have been defined with acceptable values
CALL qpolcheckpar(ok$)
IF ok$ <> "y" THEN EXIT SUB

'Calculation Of GL:
'VT2 = N * L * 1/2 * G * r^2            'N turns at r
'    = N * L * 1/2 * G * R^2 * 2        'N at R and N at -R
'    = 1/2 * N * GL * 2 * R^2
'GL = VT2 * (2 / (N * 2 * R^2))
'GL = CoilConst * VT2

'Calculation Of TH:
'-sin[n (th - TH)] = cos(n th + vtphase)
'                  = -sin(n th + vtphase - pi/2)
'-n TH = vtphase - pi/2
'TH = -(vtphase - pi/2) / n

'Get the coil integrated voltage at the second harmonic
CALL vtgetvthar2("+", vtmagpos!, vtphasepos!)
CALL vtgetvthar2("-", vtmagneg!, vtphaseneg!)
vtmag! = (vtmagpos! + vtmagneg!) / 2!
vtphase! = (vtphasepos! + vtphaseneg!) / 2!

'Compute GL
gl! = strqpolcoilconstP! * vtmag!

'Compute TH
th! = -(vtphase! - 90!) / 2!

'Log the result
CALL qpolloggl(logfileP$, strqpolcoilconstP!, vtmag!, vtphase!, gl!, th!)

END SUB

SUB qpolgetglav (gl!, sgl!, th!, sth!)

'****************************************************************************
'This subroutine has GL measured several times, then it computes the
'average and rms variation.
'
'Output:
' gl!, average of the GL measurements
' sgl!, rms variation of the GL measurements
' th!, angle of the first south pole (deg)
' sth!, rms variation of the TH measurements
'
'Zachary Wolf
'11/21/95
'****************************************************************************

'Make sure all parameters have been defined with acceptable values
CALL qpolcheckpar(ok$)
IF ok$ <> "y" THEN EXIT SUB

'Log that a strength measurement is being made
CALL qpollogstrmeas(logfileP$)

'Simplify the notation
nm% = blnmeasaveP%

'Initialize data arrays
DIM glmeas!(1 TO nm%)
DIM thmeas!(1 TO nm%)

'Message about the measurement number
PRINT
PRINT "Performing GL measurement #";

'Loop over the measurements
FOR i% = 1 TO nm%

'Message about the measurement number
PRINT USING " #"; i%;

'Measure GL
CALL qpolgetgl(gl!, th!)

'Store the results
glmeas!(i%) = gl!
thmeas!(i%) = th!

'End the loop over the measurements
NEXT i%

'End message info
PRINT

'Compute the GL average
gl! = 0!
FOR i% = 1 TO nm%
 gl! = gl! + glmeas!(i%)
NEXT i%
gl! = gl! / nm%

'Compute the GL rms deviation
sgl! = 0!
FOR i% = 1 TO nm%
 sgl! = sgl! + (glmeas!(i%) - gl!) ^ 2
NEXT i%
sgl! = SQR(sgl! / nm%)

'Compute TH average
th! = 0!
FOR i% = 1 TO nm%
 th! = th! + thmeas!(i%)
NEXT i%
th! = th! / nm%

'Compute TH rms deviation
sth! = 0!
FOR i% = 1 TO nm%
 sth! = sth! + (thmeas!(i%) - th!) ^ 2
NEXT i%
sth! = SQR(sth! / nm%)

'Print the result to the screen
PRINT "Integrated gradient GL = "; gl!; " +- "; sgl!; " T"
PRINT "South Pole Angle TH = "; th!; " +- "; sth!; " deg"

'Log the result
CALL qpollogglav(logfileP$, gl!, sgl!, th!, sth!)

END SUB

SUB qpollogcalc (logfile$)

'****************************************************************************
'This subroutine summarizes the calculations used in the program.
'
'Inputs:
' logfile$, the name of the log file
'
'Zachary Wolf
'9/22/94, 12/3/94
'****************************************************************************

'Open the log file
filenum% = FREEFILE
OPEN logfile$ FOR APPEND AS filenum%

'Header
'PRINT #filenum%,
'PRINT #filenum%, "           SUMMARY OF THE CALCULATIONS AND CONVENTIONS USED"

'Summarize the calculations used for GL vs I
'PRINT #filenum%,
'PRINT #filenum%,
'PRINT #filenum%, "Calculation Of GL:"
'PRINT #filenum%, "V2 = Nturns * velocity * B2 * L"
'PRINT #filenum%, "   = Nturns * Rcoil * ang_freq * GL * Rcoil"
'PRINT #filenum%, "   = Nturns * Rcoil^2 * ang_freq * GL"
'PRINT #filenum%, "GL = V2 / (Nturns * Rcoil^2 * ang_freq)"
'PRINT #filenum%, "For a double coil (for coil bowing correction),"
'PRINT #filenum%, "GL = ( (1 / Nturns * Rcoil^2)1 + (1 / Nturns * Rcoil^2)2 ) * (V2 / ang_freq)"
'PRINT #filenum%, "Let CoilConst = ( (1 / Nturns * Rcoil^2)1 + (1 / Nturns * Rcoil^2)2 )"
'PRINT #filenum%, "GL = CoilConst * V2 / ang_freq"

'Magnetic center
PRINT #filenum%,
PRINT #filenum%,
PRINT #filenum%, "Calculation Of The Magnetic Center:"
PRINT #filenum%, "In the quadrupole's frame, Bx' = G * y', By' = G * x'."
PRINT #filenum%, "In the coil's frame (unprimed frame) the magnetic center is at (x0, y0)."
PRINT #filenum%, "In the coil's frame,"
PRINT #filenum%, "     Bx = G * (y - y0), By = G * (x - x0)."
PRINT #filenum%, "This gives the magnetic center in terms of the measured"
PRINT #filenum%, "dipole field."
PRINT #filenum%, "     x0 = - By / G, y0 = - Bx / G"
PRINT #filenum%, "In terms of the measured integrated strengths,"
PRINT #filenum%, "     Xcenter = - (1/GL) * BL1 * sin(THspole1)"
PRINT #filenum%, "     Ycenter = - (1/GL) * BL1 * cos(THspole1)"

'Close the log file
CLOSE filenum%

END SUB

SUB qpolloggl (logfile$, coilconst!, vt2!, vtphase!, gl!, th!)

'****************************************************************************
'This subroutine writes a GL measurement to the log file.
'
'Input:
' logfile$, the name of the log file
' coilconst!, coil constant
' vt2!, coil integrated voltage at the second harmonic
' vtphase!, coil integrated voltage phase angle of the second harmonic (deg)
' gl!, integrated gradient (T)
' th!, south pole angle (deg)
'
'Zachary Wolf
'11/19/95, 10/23/96
'****************************************************************************

'Open the log file
filenum% = FREEFILE
OPEN logfile$ FOR APPEND AS filenum%

'Print the results
PRINT #filenum%, TIME$; "  Integrated Gradient Measurement:"
PRINT #filenum%, "          CoilConst = "; coilconst!; " 1/m^2"
PRINT #filenum%, "          VT2 = "; vt2!; " Vs, VTphase = "; vtphase!; " deg"
PRINT #filenum%, "          GL = "; gl!; " T, TH = "; th!; " deg"

'Close the log file
CLOSE filenum%

END SUB

SUB qpollogglav (logfile$, gl!, sgl!, th!, sth!)

'****************************************************************************
'This subroutine writes a summary of an average GL measurement to the log file.
'
'Input:
' logfile$, the name of the log file
' gl!, the average of the integrated gradient measurements
' sgl!, the rms variation of gl
' th!, the average of the south pole angle measurements (deg)
' sth!, the rms variation of TH
'
'Zachary Wolf
'11/19/95
'****************************************************************************

'Open the log file
filenum% = FREEFILE
OPEN logfile$ FOR APPEND AS filenum%

'Print the results
PRINT #filenum%, TIME$; "  Average Integrated Gradient:"
PRINT #filenum%, "          GL = "; gl!; " +- "; sgl!; " T"
PRINT #filenum%, "          TH = "; th!; " +- "; sth!; " deg"

'Close the log file
CLOSE filenum%

END SUB

SUB qpollogstrmeas (logfile$)

'****************************************************************************
'This subroutine records a measurement to the log file.
'
'Input:
' logfile$, the name of the log file
'
'Zachary Wolf
'11/25/95
'****************************************************************************

'Open the log file
filenum% = FREEFILE
OPEN logfile$ FOR APPEND AS filenum%

'Print the results
PRINT #filenum%,
PRINT #filenum%, TIME$; "  Magnet Strength Measurement..."

'Close the log file
CLOSE filenum%

END SUB

SUB qpolmaingl

'****************************************************************************
'This subroutine controls the magnetic measurements of a quadrupole.
'
'Zachary Wolf
'1/4/95, 5/6/95, 8/19/96
'****************************************************************************

'Make sure all parameters have been defined with acceptable values
CALL qpolcheckpar(ok$)
IF ok$ <> "y" THEN EXIT SUB

'Measure all temperatures
CALL pmtgettslice(name1$, t1!, name2$, t2!, name3$, t3!)
CALL pmtdattslice(strdatfileP$, name1$, t1!, name2$, t2!, name3$, t3!)
CALL pmtdattslice(hardatfileP$, name1$, t1!, name2$, t2!, name3$, t3!)

'Message
PRINT
PRINT "Beginning the quadrupole field measurement cycle..."

'Simplify the notation
nh% = nhardisplayP%

'Initialize data arrays
DIM bln!(1 TO nh%)                'field strength at Rcoil of the n'th harmonic
DIM sbln!(1 TO nh%)               'standard deviation of Bn
DIM thn!(1 TO nh%)                'field angle
DIM sthn!(1 TO nh%)               'standard deviation of THn
DIM rn!(1 TO nh%)                 'BLn/BL2
DIM srn!(1 TO nh%)                'standard deviation of rn

'Dipole Strength (for IR magnets)

'Message
PRINT
PRINT "Measuring the dipole strength..."

'Connect the dipole coil to the integrator
CALL muxsetchan("dipole")

'Set the integrator gain
CALL vtautogain

'Get the integrated field strength
CALL dpolgetblav(bl!, sbl!, th!, sth!)

'Write the strength measurement to the data file
CALL dpoldatbl(strdatfileP$, bl!, sbl!, th!, sth!)

'Write the strength measurement to the plot file
CALL dpolpltbl(strpltfileP$, bl!, sbl!, th!, sth!)

'Quadrupole Strength

'Message
PRINT
PRINT "Measuring the quadrupole strength..."

'Connect the quadrupole coil to the integrator
CALL muxsetchan("quadrupole")

'Set the integrator gain
CALL vtautogain

'Get the integrated gradient
CALL qpolgetglav(gl!, sgl!, th!, sth!)

'Write the strength measurement to the data file
CALL qpoldatgl(strdatfileP$, gl!, sgl!, th!, sth!)

'Write the strength measurement to the plot file
CALL qpolpltgl(strpltfileP$, gl!, sgl!, th!, sth!)

'Harmonics

'Message
PRINT
PRINT "Measuring BL1 and BL2 with the main winding..."

'Connect the main coil to the integrator
CALL muxsetchan("main")

'Set the integrator gain
CALL vtautogain

'Get the field harmonics at the coil radius
CALL blgetblnav("main", nh%, bln!(), sbln!(), thn!(), sthn!())

'Keep the dipole and quadrupole terms
bl1! = bln!(1)
sbl1! = sbln!(1)
th1! = thn!(1)
sth1! = sthn!(1)
bl2! = bln!(2)
sbl2! = sbln!(2)
th2! = thn!(2)
sth2! = sthn!(2)

'Message
PRINT
PRINT "Measuring BLn with the bucking winding..."

'Connect the main + dipole + quadrupole bucking coils to the integrator
CALL muxsetchan("dqbuck")

'Set the integrator gain
CALL vtautogain

'Get the field harmonics at the coil radius
CALL blgetblnav("dqbuck", nh%, bln!(), sbln!(), thn!(), sthn!())

'Fill in the dipole and quadrupole terms from using the main winding
bln!(1) = bl1!
sbln!(1) = sbl1!
thn!(1) = th1!
sthn!(1) = sth1!
bln!(2) = bl2!
sbln!(2) = sbl2!
thn!(2) = th2!
sthn!(2) = sth2!

'Compute the strength ratios of the harmonics
CALL blcalcratio(nh%, bln!(), sbln!(), nmain%, rn!(), srn!())

'Compute the location of the magnetic center
CALL qpolcalcxy(nh%, bln!(), sbln!(), thn!(), sthn!(), x!, sx!, y!, sy!)

'Write the harmonics to the data file
CALL bldatblhar(hardatfileP$, nmain%, nh%, bln!(), sbln!(), thn!(), sthn!(), rn!(), srn!())

'Write the magnetic center to the data file
CALL qpoldatxy(hardatfileP$, thn!(), x!, sx!, y!, sy!)

'Write the harmonics to the plot file
CALL blpltblhar(harpltfileP$, nh%, rn!(), srn!())

'Measure all temperatures again, after the coil measurements
CALL pmtgettslice(name1$, t1!, name2$, t2!, name3$, t3!)
CALL pmtdattslice(strdatfileP$, name1$, t1!, name2$, t2!, name3$, t3!)
CALL pmtdattslice(hardatfileP$, name1$, t1!, name2$, t2!, name3$, t3!)

'Write a summary of the calculations to the GL data file
CALL bllogcalc(strdatfileP$)
CALL qpollogcalc(strdatfileP$)

'Write a summary of the calculations to the harmonics data file
CALL bllogcalc(hardatfileP$)
CALL qpollogcalc(hardatfileP$)

END SUB

SUB qpolmainglvsi

'****************************************************************************
'This subroutine controls the magnetic measurements of a quadrupole at
'various currents.
'
'Zachary Wolf
'1/4/95, 5/6/95
'****************************************************************************

'Make sure all parameters have been defined with acceptable values
CALL qpolcheckpar(ok$)
IF ok$ <> "y" THEN EXIT SUB

'Message
PRINT
PRINT "Beginning the quadrupole field measurement cycle..."

'Simplify the notation
ni% = nimagtestP%
nh% = nhardisplayP%

'Initialize data arrays
DIM bln!(1 TO nh%)                'field strength at Rcoil of the n'th harmonic
DIM sbln!(1 TO nh%)               'standard deviation of Bn
DIM thn!(1 TO nh%)                'field angle
DIM sthn!(1 TO nh%)               'standard deviation of THn
DIM rn!(1 TO nh%)                 'BLn/BL2
DIM srn!(1 TO nh%)                'standard deviation of rn

'Loop over test currents
FOR i% = 1 TO ni%

'Ramp to the desired current
CALL imagramp(imagtestP!(i%))

'Measure the current
CALL imaggetiav(imag!, simag!)

'Connect the quadrupole coil to the integrator
CALL muxsetchan("quadrupole")

'Set the integrator gain
CALL vtautogain

'Get the integrated gradient
CALL qpolgetglav(gl!, sgl!, th!, sth!)

'Compute the transfer function
tf! = gl! / imag!
stf! = tf! * SQR((sgl! / gl!) ^ 2 + (simag! / imag!) ^ 2)

'Write the strength measurement to the data file
CALL qpoldatglvsi(strdatfileP$, imag!, simag!, gl!, sgl!, tf!, stf!)

'Write the strength measurement to the plot file
CALL qpolpltglvsi(strpltfileP$, imag!, gl!, sgl!, tf!, stf!)

'Do the following only if the harmonics are desired
IF imagharflagP!(i%) > 0 THEN

 'Connect the main coil to the integrator
 CALL muxsetchan("main")

 'Set the integrator gain
 CALL vtautogain

 'Get the field harmonics at the coil radius
 CALL blgetblnav("main", nh%, bln!(), sbln!(), thn!(), sthn!())

 'Compute the strength ratios of the harmonics
 CALL blcalcratio(nh%, bln!(), sbln!(), nmain%, rn!(), srn!())

 'Compute the location of the magnetic center
 CALL qpolcalcxy(nh%, bln!(), sbln!(), thn!(), sthn!(), x!, sx!, y!, sy!)

 'Write the harmonics to the data file
 CALL bldatblvsihar(hardatfileP$, imag!, simag!, nmain%, nh%, bln!(), sbln!(), thn!(), sthn!(), rn!(), srn!())

 'Write the magnetic center to the data file
 CALL qpoldatxy(hardatfileP$, thn!(), x!, sx!, y!, sy!)

 'Write the harmonics to the plot file
 legend$ = "I = " + STR$(imagtestP!(i%)) + " A"
 CALL blpltblvsihar(harpltfileP$, legend$, nh%, rn!(), srn!())

'End harmonics check
END IF

'End loop over currents
NEXT i%

'Write a summary of the calculations to the GL data file
CALL filenewpage(strdatfileP$)
CALL qpollogcalc(strdatfileP$)

'Write a summary of the calculations to the harmonics data file
CALL filenewpage(hardatfileP$)
CALL qpollogcalc(hardatfileP$)

END SUB

SUB qpolpltgl (pltfile$, gl!, sgl!, th!, sth!)

'****************************************************************************
'This subroutine write a summary of the GL measurements for plotting.
'
'Input:
' pltfile$, the name of the plot file
' gl!, the integrated gradient at the measured current
' sgl!, the standard deviation of gl
' th!, south pole angle
' sth!, rms variation of TH
'
'Zachary Wolf
'11/19/95, 8/19/96
'****************************************************************************

'Keep track of Ncall for the header
STATIC ncall%
ncall% = ncall% + 1

'Open the plot file
filenum% = FREEFILE
OPEN pltfile$ FOR APPEND AS filenum%

'Write the title for plotting on the first call
IF ncall% = 1 THEN
 PRINT #filenum%, "/et g '" + magnameP$ + ", Run " + runP$ + "'"
 PRINT #filenum%, ";gl sgl th sth"
END IF

'Write the values
PRINT #filenum%, USING " ##.#####"; gl!;
PRINT #filenum%, USING " ##.#####"; sgl!;
PRINT #filenum%, USING " ##.#####"; th!;
PRINT #filenum%, USING " ##.#####"; sth!

'Close the log file
CLOSE filenum%

END SUB

SUB qpolpltglvsi (pltfile$, imag!, gl!, sgl!, tf!, stf!)

'****************************************************************************
'This subroutine write a summary of the GL vs I measurements for plotting.
'
'Input:
' pltfile$, the name of the plot file
' imag!, the measured current
' gl!, the integrated gradient at the measured current
' sgl!, the standard deviation of gl
' tf!, the transfer function at the measured current
' stf!, the standard deviation of tf
'
'Zachary Wolf
'11/19/95
'****************************************************************************

'Keep track of Ncall for the header
STATIC ncall%
ncall% = ncall% + 1

'Open the plot file
filenum% = FREEFILE
OPEN pltfile$ FOR APPEND AS filenum%

'Write the title for plotting on the first call
IF ncall% = 1 THEN
 PRINT #filenum%, "/et g '" + magnameP$ + ", Run " + runP$ + "'"
 PRINT #filenum%, ";i gl sgl tf stf"
END IF

'Write the values
PRINT #filenum%, USING "####.#"; imag!;
PRINT #filenum%, USING " ##.#####"; gl!;
PRINT #filenum%, USING " ##.#####"; sgl!;
PRINT #filenum%, USING " ###.#####"; tf!;
PRINT #filenum%, USING " ###.#####"; stf!

'Close the log file
CLOSE filenum%

END SUB

SUB qpolsetpar (logfile$, strdatfile$, strpltfile$, hardatfile$, harpltfile$, magname$, run$, strqpolcoilconst!, blnmeasave%, nimagtest%, imagtest!(), imagharflag!(), nhardisplay%)

'****************************************************************************
'This subroutine sets the parameters required by this module.
'
'Input:
' logfile$, path and name of the log file
' strdatfile$, magnet strength data file
' strpltfile$, magnet strength plot file
' hardatfile$, harmonics data file
' harpltfile$, harmonics plot file
' magname$, magnet name
' run$, run number
' strqpolcoilconst!, strength coil quadrupole coil constant
' blnmeasave%, number of measurements for averaging
' nimagtest%, number of test currents
' imagtest!(), test currents
' imagharflag!(), flag for harmonics measurement
' nhardisplay%, number of harmonics to record
'
'Zachary Wolf
'11/20/95
'****************************************************************************

'Set the parameters
logfileP$ = logfile$
strdatfileP$ = strdatfile$
strpltfileP$ = strpltfile$
hardatfileP$ = hardatfile$
harpltfileP$ = harpltfile$
magnameP$ = magname$
runP$ = run$
strqpolcoilconstP! = strqpolcoilconst!
blnmeasaveP% = blnmeasave%
nimagtestP% = nimagtest%
DIM imagtestP!(1 TO nimagtestP%)
DIM imagharflagP!(1 TO nimagtestP%)
FOR i% = 1 TO nimagtestP%
 imagtestP!(i%) = imagtest!(i%)
 imagharflagP!(i%) = imagharflag!(i%)
NEXT i%
nhardisplayP% = nhardisplay%

END SUB