; NAME:
;     originally APPLY_CAL, renamed calib_apply
; PURPOSE:
;     Applies the total power, or amplitude and phase calibration contained
;     in REFCAL, to a dynamic spectrum for a single antenna or baseline and
;     a single polarization.
; CATEGORY:
;     OVRO APC CALIBRATION
; CALLING SEQUENCE:
;     apply_cal,amp,pha,time,freq,ipol,ibl,refcal
; INPUTS:
;     amp       a dynamic spectrum of size (NF,NTIMES), for a single antenna
;                 (in the case of total power), or a single baseline (in the
;                 case of amplitude data), and a single polarization.
;                 NB: the updated values are returned in AMP, overwriting
;                 its contents.
;     pha       an array of the same type as AMP, containing the phase for
;                 the same baseline.  If AMP is total power (as indicated by
;                 IBL value) then pha is not used, but must be supplied (e.g.
;                 use a constant such as zero).
;                 NB: the updated values are returned in PHA, overwriting
;                 its contents.
;     time      an array of size (NTIMES) containing the times for the
;                 AMP and PHA arrays [msec].
;     freq      an array of size (NF) containing the frequencies for the
;                 AMP and PHA arrays [GHz].
;     ipol      a constant specifying the polarization type: 0=RCP, 1=LCP,
;                 and 2=LIN, of the data in AMP and PHA.
;     ibl       a constant specifying the antenna or baseline of the data in
;                 AMP and PHA.  Note that this is the index into the REFCAL
;                 arrays, and its meaning will depend on the number of antennas
;                 and baselines.  Values 0 to NANT-1 indicate total power for
;                 antenna ANT[IBL].  Values NANT to NANT+NBL-1 indicate
;                 amplitude and phase for baseline IBL-NANT.  For example, here
;                 are the 21 values for 6 antennas: 0=Ant 1, 1=Ant 2, 2=Ant 4,
;                 3=Ant 5, 4=Ant 6, 5=Ant 7, 6=Baseline 12, 7=Baseline 14,
;                 8=Baseline 15, ... 11=Baseline 24, 12=Baseline 25, ...
;                 14=Baseline 27, 15=Baseline 45, ... 20=Baseline 67.
;     refcal    the complete REFCAL structure as returned from NEWSCAN.  The
;                 information from REFCAL is applied to the AMP and PHA arrays
;                 to calibrate the total power or amplitudes and phases.  In
;                 the case of phase, both the reference calibration and the
;                 daily (time-dependent) phase calibration are applied.  In the
;                 case of total power or amplitude, both the reference calibration
;                 and the update calibration are applied.
; OPTIONAL (KEYWORD) INPUT PARAMETERS:
; ROUTINES CALLED:
;     lobe
; OUTPUTS:
;     amp       the same dynamic spectrum as the input, except with calibration
;                 applied.  NB: this overwrites the contents of the input.
;     pha       the same dynamic spectrum as the input, except with calibration
;                 applied.  NB: this overwrites the contents of the input.
; COMMENTS:
; SIDE EFFECTS:
; RESTRICTIONS:
; MODIFICATION HISTORY:
;     Written 29-Jul-2000 by Dale E. Gary
;     07-Oct-2000  DG
;       Fixed problem with application of phase calibration for a
;       single frequency.
;-

pro calib_apply,amp,pha,time,freq,ipol,ibl,refcal

   ; Compare input frequency list with frequencies in REFCAL
   nfreq = n_elements(freq)
   ntimes = n_elements(time)
   nant = refcal.nant

   ; This should work for one, or many, specified frequencies.
   ; IFREQ will be an array of indexes into the frequency array.
   ifreq = intarr(nfreq)
   for i = 0, nfreq-1 do begin
      ; Find indexes into frequency list for all frequencies requested.
      ; These should match exactly, but since they are REAL, compare
      ; them scaled by factor of 10
      ifreq[i] = where(nint(freq[i]*10) eq nint((*refcal.pfrq)*10))
   endfor

   ; If there were any mismatches, bail out
   bad = where(ifreq eq -1,nbad)
   if (nbad ne 0) then begin
      print,'APPLY_CAL: Invalid frequency list.  No calibration applied.'
      return
   endif

   ; Get amplitude factors (independent of time) for either baseline or
   ; total power.

   ; Update factors can be missing, in which case they are NaN. However,
   ; a missing update factor ought to be 1.0, so set them that way...
   updfac = (*refcal.pfacupd)[ipol,ibl,ifreq]
   bad = where(finite(updfac) eq 0,nbad)
   if (nbad ne 0) then updfac[bad] = 1.0

   ampfac = (*refcal.pfactors)[ipol,ibl,ifreq]*updfac
   ; Apply it
   amp = amp/reform(reform(ampfac)#replicate(1,ntimes))

   ; In the case of baseline data, we also need the (time dependent) phase
   ; corrections.
   if (ibl ge nant) then begin
      ; Reference calibration part (time independent)
      phase = reform((*refcal.pphase)[ipol,ibl-nant,ifreq])

      ; Time dependent part (from DAILY calibration)
      t1 = refcal.phaseupd[0].tref
      t2 = refcal.phaseupd[1].tref
      po1 = (*refcal.phaseupd[0].ppoff)[ibl-nant]
      po2 = (*refcal.phaseupd[1].ppoff)[ibl-nant]
      ps1 = (*refcal.phaseupd[0].pdpdf)[ibl-nant]
      ps2 = (*refcal.phaseupd[1].pdpdf)[ibl-nant]
      if (t1 eq t2) then tfac = 0 else tfac = double(time-t1)/(t2-t1)

      ; Apply time-dependent part.  Result is of size (NTIMES,NFREQ)
      if (nfreq eq 1) then begin
         phase = phase#replicate(1,ntimes) + (po1 + tfac*lobe(po2-po1,/mid0)) + freq*(ps1 + tfac*(ps2-ps1))
      endif else begin
         phase = phase#replicate(1,ntimes) $
               + replicate(1,nfreq)#(po1 + tfac*lobe(po2-po1,/mid0)) $
               + freq#(ps1 + tfac*(ps2-ps1))
      endelse

      ; Apply phase corrections
      pha = lobe(pha - phase)

   endif

return
end