;+
; Combine (calibrate) a frequency switched integration.  The resulting
; data container is returned.
;
; @param sig {in}{required}{type=spectrum} The signal spectrum with no
; cal.
; @param ref {in}{required}{type=spectrum} The reference spectrum with
; no cal.
; @param sigwcal {in}{required}{type=spectrum} The signal with cal on spectrum.
; @param refwcal {in}{required}{type=spectrum} The reference with cal on spectrum.
;
; @returns A spectrum data container containing the result.
;
; @examples
; <pre>
; get, mc_scan=22, cal='F', sig='T', pol='XX', if_num=1, int=1 ; sig
; copy,0,10 
; get, mc_scan=22, cal='T', sig='T', pol='XX', if_num=1, int=1 ;sigwcal
; copy,0,11
; get, mc_scan=22, cal='F', sig='F', pol='XX', if_num=1, int=1 ; ref
; copy,0,12
; get, mc_scan=22, cal='T', sig='F', pol='XX', if_num=1, int=1 ; refwcal
; copy,0,13
; result = fscal(!g.s[10],!g.s[11], !g.s[12], !g.s[13])
; set_data_container, result, 0 ; copy result back to buffer 0 and show the result
; data_free, result  ; clean up pointers if result is no longer necessary
; </pre>
;
; @uses <a href="data_valid.html">data_valid</a>
; @uses <a href="data_valid.html">data_copy</a>
;
; @version $Id: fscal.pro,v 1.2 2005/02/28 18:28:39 jbraatz Exp $
;-
function fscal, sig, ref, sigwcal, refwcal
    compile_opt idl2

    ; are the arguments valid and have the same length
    nsig = data_valid(sig)
    nref = data_valid(ref)
    nsigwcal = data_valid(sigwcal)
    nrefwcal = data_valid(refwcal)

    if (nsig le 0 or nref le 0 or nsigwcal le 0 or nrefwcal le 0) then begin
        message, 'One or more of the data containers are invalid'
    endif

    if (nsig ne nref or nsig ne nsigwcal or nsig ne nrefwcal) then begin
        message, 'One or more of the data containers does not have the same length'
    endif

    nchans = n_elements(*sig.data_ptr)
    pct10 = nchans/10
    pct90 = nchans - pct10
    tsysSig = (mean((*sigwcal.data_ptr)[pct10:pct90] + (*sig.data_ptr)[pct10:pct90]) / $
               mean((*sigwcal.data_ptr)[pct10:pct90] - (*sig.data_ptr)[pct10:pct90]))*sig.mean_tcal - $
               sig.mean_tcal/2.0
    sigResult = (*sigwcal.data_ptr + *sig.data_ptr)/2.0

    tsysRef = (mean((*refwcal.data_ptr)[pct10:pct90] + (*ref.data_ptr)[pct10:pct90]) / $
               mean((*refwcal.data_ptr)[pct10:pct90] - (*ref.data_ptr)[pct10:pct90]))*ref.mean_tcal - $
               ref.mean_tcal/2.0
    refResult = (*refwcal.data_ptr + *ref.data_ptr)/2.0

    data_copy, sig, result

    *result.data_ptr = ((sigResult - refResult) / refResult) * tsysRef
    result.tsys = tsysRef
    result.tsysref = tsysSig
    result.freq_switch_offset = ref.reference_frequency - sig.reference_frequency + $
                sig.frequency_interval * (sig.reference_channel - ref.reference_channel)

    return, result
end