DECLARE SUB qwcalcgl (nx%, x0!(), dx!, vtvsx!(), svtvsx!(), gl!, sgl!)
DECLARE SUB qwcheckpar (ok$)
DECLARE SUB qwdatglvsi (logfile$, imag!, simag!, gl!, sgl!)
DECLARE SUB qwgetgl (gl!, sgl!)
DECLARE SUB qwlinfit (n%, x!(), y!(), m!, sm!, b!, sb!)
DECLARE SUB qwloggl (logfile$, gl!, sgl!)
DECLARE SUB qwpltglvsi (logfile$, imag!, simag!, gl!, sgl!)
DECLARE SUB vtwgetvtvsx (nx%, x0!(), dx!, vtvsx!(), svtvsx!())
DECLARE SUB vtwmoveabs (x!)
DECLARE SUB imaggetiav (imag!, simag!)
DECLARE SUB imagramp (inom!)

'****************************************************************************
'Module QW.BAS
'
'These subroutines perform calculations related to quadrupole magnet
'stretched wire measurements.
'
'Zachary Wolf
'1/14/95
'****************************************************************************

'Common area for shared parameters
COMMON SHARED /qwire/ logfileP$, datfileP$, pltfileP$, magnameP$, runP$, dxcmP!, nx0P%, x0P!(), nturnsP%, nimagtestP%, imagtestP!()

SUB qwcalcgl (nx%, x0!(), dx!, vtvsx!(), svtvsx!(), gl!, sgl!)

'****************************************************************************
'This subroutine calculates the integrated gradient from the measurements
'of the integrated voltage as a function of x.
'
'Input:
' nx%, the number of x positions
' x0!(1 to nx%), the position of the center of each wire motion
' dx!, the distance the wire moves for each motion
' vtvsx!(1 to nx%), the integrated voltage measurements
' svtvsx!(1 to nx%), the standard deviation on each vt
'
'Output:
' gl!, the integrated gradient of the magnet
' sgl!, the standard deviation on gl
'
'Zachary Wolf
'1/14/95
'****************************************************************************

'Fit the vt vs x data points with a line
CALL qwlinfit(nx%, x0!(), vtvsx!(), m!, sm!, b!, sb!)   'm=dVT/dX0  (Vs/cm)

'Compute the integrated gradient
'gl = (1/delX)*(dVT/dX0)
m! = m! * 100!                                '(Vs/m)
m! = m! / nturnsP%                            'take out Nturns, VT now flux change in magnet
gl! = m! / dx!                                'gl = (1/delx)(dVT/dX0), ((1/cm)(Tm^2/m))
gl! = gl! * 100!                              '(T)

'Error
sm! = sm! * 100!                              '(Vs/m)
sm! = sm! / nturnsP%                          'take out Nturns, VT now flux change in magnet
sgl! = sm! / dx!                              'gl = (1/delx)(dVT/dX0), ((1/cm)(Tm^2/m))
sgl! = sgl! * 100!                            '(T)

END SUB

SUB qwcheckpar (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"

'logfileP$
IF logfileP$ = "" THEN
 PRINT "QWIRE: log file not defined"
 EXIT SUB
END IF

'datfileP$
IF logfileP$ = "" THEN
 PRINT "QWIRE: dat file not defined"
 EXIT SUB
END IF

'pltfileP$
IF logfileP$ = "" THEN
 PRINT "QWIRE: plt file not defined"
 EXIT SUB
END IF

'magnameP$
IF magnameP$ = "" THEN
 PRINT "QWIRE: magnet name not defined"
 EXIT SUB
END IF

'runP$
IF logfileP$ = "" THEN
 PRINT "QWIRE: run number not defined"
 EXIT SUB
END IF

'dxcmP!
IF dxcmP! <= 0! OR dxcmP! > 10! THEN
 PRINT "QWIRE: dxcm has improper value"
 EXIT SUB
END IF

'nx0P%
IF nx0P% <= 0 OR nx0P% > 50 THEN
 PRINT "QWIRE: nx0 has improper value"
 EXIT SUB
END IF

'nturnsP%
IF nturnsP% <= 0 OR nturnsP% > 200 THEN
 PRINT "QWIRE: nturns has improper value"
 EXIT SUB
END IF

'nimagtestP%
IF nimagtestP% <= 0 OR nimagtestP% > 50 THEN
 PRINT "QWIRE: nimagtest has improper value"
 EXIT SUB
END IF

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

END SUB

SUB qwdatglvsi (logfile$, imag!, simag!, gl!, sgl!)

'****************************************************************************
'This subroutine writes the calculated integrated gradient of the magnet.
'
'Input:
' logfile$, the name of the log file
' gl!, the integrated gradient of the magnet
' sgl!, the standard deviation of GL
'
'Zachary Wolf
'1/14/95
'****************************************************************************

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

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

'Write the header
IF ncall% = 1 THEN
 PRINT #filenum%,
 PRINT #filenum%, "       INTEGRATED GRADIENT VS CURRENT"
 PRINT #filenum%,
 PRINT #filenum%, "  Imag       sImag        GL         sGL"
 PRINT #filenum%, "   (A)        (A)         (T)        (T)"
 PRINT #filenum%, "---------+----------  ----------+-----------"
END IF

'Calculated Results
PRINT #filenum%, USING "#####.###"; imag!;
PRINT #filenum%, USING "  #####.###"; simag!;
PRINT #filenum%, USING "  ####.#####"; gl!;
PRINT #filenum%, USING "  ####.#####"; sgl!

'Close the log file
CLOSE filenum%

END SUB

SUB qwgetgl (gl!, sgl!)

'****************************************************************************
'This subroutine supervises the wire measurements and integrate gradient
'calculation.
'
'Output:
' gl!, the integrated gradient
' sgl!, the standard deviation of gl!
'
'Zachary Wolf
'10/31/94
'****************************************************************************

'Initialize data arrays
DIM vtvsx!(1 TO nx0P%)              'average voltage integral at each position
DIM svtvsx!(1 TO nx0P%)             'standard deviation of the voltage integral measurements

'Message
PRINT
PRINT "Beginning the stretched wire measurement cycle..."

'Measure the voltage integral at points across the magnet aperture
CALL vtwgetvtvsx(nx0P%, x0P!(), dxcmP!, vtvsx!(), svtvsx!())

'Plot VT vs X
'CALL qwpltvtvsx(pltfileP$, magnameP$, runP$, nx0P%, x0P!(), vtvsx!(), svtvsx!())

'Calculate the integradient of the quadrupole
CALL qwcalcgl(nx0P%, x0P!(), dxcmP!, vtvsx!(), svtvsx!(), gl!, sgl!)

'Message
PRINT
PRINT "The integrated gradient is GL = "; gl!; " +- "; sgl!; " T"

'Write the result to the log file
CALL qwloggl(logfileP$, gl!, sgl!)

'Message
PRINT
PRINT "Moving back to the home position..."

'Move the wire back to the home position
CALL vtwmoveabs(0!)

END SUB

SUB qwglvsi

'****************************************************************************
'This subroutine loops over the test currents and performs a wire scan
'at each current.
'
'Zachary Wolf
'5/16/95
'****************************************************************************

'Make sure all parameters have been initialized
CALL qwcheckpar(ok$)
IF ok$ <> "y" THEN
 PRINT "QW: the parameters in qwblvsi have a problem"
 EXIT SUB
END IF

'Loop over currents
FOR i% = 1 TO nimagtestP%

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

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

'Perform the stretched wire measurement
CALL qwgetgl(gl!, sgl!)

'Record the result in the data file
CALL qwdatglvsi(datfileP$, imag!, simag!, gl!, sgl!)

'Plot the results
CALL qwpltglvsi(pltfileP$, imag!, simag!, gl!, sgl!)

'End current loop
NEXT i%

END SUB

SUB qwlinfit (n%, x!(), y!(), m!, sm!, b!, sb!)

'****************************************************************************
'This subroutine does a least squares linear fit to the (x,y) data points.
'It returns the slope and intercept of the line with errors.
'See Bevington, pp. 104, 114
'
'Input:
' n%, the number of (x,y) points
' x!(1 to n%), the x values
' y!(1 to n%), the y values
'
'Output:
' m!, the fitted line's slope
' sm!, error estimate for m
' b!, the fitted line's intercept
' sb!, error estimate for b
'
'Zachary Wolf
'5/28/95
'****************************************************************************

'Compute various sums
SUMxy! = 0!
SUMxx! = 0!
SUMx! = 0!
SUMy! = 0!
FOR i% = 1 TO n%
 SUMxy! = SUMxy! + x!(i%) * y!(i%)
 SUMxx! = SUMxx! + x!(i%) * x!(i%)
 SUMx! = SUMx! + x!(i%)
 SUMy! = SUMy! + y!(i%)
NEXT i%

'Compute the determinant
d! = n% * SUMxx! - SUMx! * SUMx!

'Compute the slope and intercept
m! = (n% * SUMxy! - SUMx! * SUMy!) / d!
b! = (SUMxx! * SUMy! - SUMxy! * SUMx!) / d!

'Compute the errors

'Estimate the error on each y value
sigsq! = 0!
FOR i% = 1 TO n%
 sigsq! = sigsq! + (y!(i%) - (m! * x!(i%) + b!)) ^ 2
NEXT i%
sigsq! = sigsq! / (n% - 2)

'Use the estimated y error to find sb! and sm!
sb! = SQR(sigsq! * SUMxx / d!)
sm! = SQR(n% * sigsq! / d!)

END SUB

SUB qwloggl (logfile$, gl!, sgl!)

'****************************************************************************
'This subroutine writes a GL measurement result to the log file.
'
'Input:
' logfile$, the name of the log file
' gl!, the integrated gradient of the magnet
' sgl!, the standard deviation of GL
'
'Zachary Wolf
'12/18/95
'****************************************************************************

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

'Write the result
PRINT #filenum%, TIME$; "  The integrated gradient is GL = "; gl!; " +- "; sgl!; " T"

'Close the log file
CLOSE filenum%

END SUB

SUB qwpltglvsi (logfile$, imag!, simag!, gl!, sgl!)

'****************************************************************************
'This subroutine writes the calculated integrated gradient of the magnet.
'
'Input:
' logfile$, the name of the log file
' gl!, the integrated gradient of the magnet
' sgl!, the standard deviation of GL
'
'Zachary Wolf
'1/14/95
'****************************************************************************

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

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

'Write the header
IF ncall% = 1 THEN
 PRINT #filenum%, ";I sI GL sGL"
END IF

'Calculated Results
PRINT #filenum%, USING "#####.###"; imag!;
PRINT #filenum%, USING "  #####.###"; simag!;
PRINT #filenum%, USING "  ##.######"; gl!;
PRINT #filenum%, USING "  ##.######"; sgl!

'Close the log file
CLOSE filenum%

END SUB

SUB qwpltvtvsx (logfile$, magname$, run$, nx%, x!(), vt!(), svt!())

'****************************************************************************
'This subroutine write a summary of the VT vs X measurements for plotting.
'
'Input:
' logfile$, the name of the log file
' magname$, the magnet name
' run$, the run number
' nx%, the number of x positions
' x!(1 to nx%), the x positions
' vt!(1 to nx%), the integrated voltage at each x position
' svt!(1 to nx%), the standard deviation of vt
'
'Zachary Wolf
'12/27/94
'****************************************************************************

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

'Write the title for plotting
PRINT #filenum%, "/et g '" + magname$ + ", Run " + run$ + "'"
PRINT #filenum%, ";x vt svt"

'Write the values
FOR i% = 1 TO nx%
 PRINT #filenum%, USING "##.####"; x!(i%);
 PRINT #filenum%, USING " ###.#######"; vt!(i%);
 PRINT #filenum%, USING " ###.#######"; svt!(i%)
NEXT i%

'Close the log file
CLOSE filenum%

END SUB

SUB qwsetpar (logfile$, datfile$, pltfile$, magname$, run$, dxcm!, nx0%, x0!(), nturns%, nimagtest%, imagtest!())

'****************************************************************************
'This subroutine sets parameters used in this module.
'
'Input:
' logfile$, name of the log file
' datfile$, name of the data file
' pltfile$, name of the plot file
' magname$, magnet name
' run$, run number
' dxcm!, distance to move the wire in a scan
' nx0%, number of x-positions at which to do wire scans
' x0!(1 to nx0%), array of x-postions
' nturnsP%, number of turns in the wire bundle
' nimagtest%, number of test currents
' imagtest!(1 to nimagtest%), test currents
'
'Zachary Wolf
'9/21/95
'****************************************************************************

'Save the parameters in the local common block for future use
logfileP$ = logfile$
datfileP$ = datfile$
pltfileP$ = pltfile$
magnameP$ = magname$
runP$ = run$
dxcmP! = dxcm!
nx0P% = nx0%
DIM x0P!(1 TO nx0P%)
FOR i% = 1 TO nx0P%
 x0P!(i%) = x0!(i%)
NEXT i%
nturnsP% = nturns%
nimagtestP% = nimagtest%
DIM imagtestP!(1 TO nimagtestP%)
FOR i% = 1 TO nimagtestP%
 imagtestP!(i%) = imagtest!(i%)
NEXT i%

END SUB