PRO plot_t_fun, ndata, trial, t_fun, plot_title, plot_file_path, plot_filename
; Plots a transfer function with ndata channels taken from
; data of trial time duration
;
; Author:  Robert Kessel
;          Electro-Optics Technology Section
;          Code 8123
;          Naval Research Laboratory
;          Washington, D.C.  20375-5354
; email:   kessel@ncst.nrl.navy.mil
;
; This code was written for Bill Palya's, Don Walter's, and Nagendra
; Tripathi's use at JSU.  However, anyone else who finds it useful
; should feel free to use it, rewrite it, or whatever else, save
; charging any money, goods or services from any other person for
; the code or use thereof.
; (5/22/00)
;

; Create some arrays
half_freqs=FINDGEN((ndata/2) + 1)
freqs=FLTARR(ndata)

; Make the frequency scale
half_freqs = half_freqs/trial
freqs[0:(ndata/2)] = half_freqs[0:(ndata/2)]
FOR i=((ndata/2)+1),(ndata-1) DO BEGIN
  freqs[i] = -half_freqs[ndata-i]
ENDFOR

; Unwrap FFT array and generate phase array
ndata = 256
trial = 2000.
unwrap_fft_freq, t_fun, ndata, trial, freqs, freqs_v2, fft_v2
milli_freqs_v2 = 1000.*freqs_v2
phase = FLTARR(256)
FOR i=0,255 DO BEGIN
  IF ((FLOAT(fft_v2[i]) EQ 0.) AND $
        (IMAGINARY(fft_v2[i]) GT 0.)) THEN BEGIN
    phase[i] = !PI/2.
  ENDIF ELSE BEGIN
    IF ((FLOAT(fft_v2[i]) EQ 0.) AND $
      (IMAGINARY(fft_v2[i]) LT 0.)) THEN BEGIN
    phase[i] = -!PI/2.
    ENDIF ELSE BEGIN
      IF ((FLOAT(fft_v2[i]) EQ 0.) AND $
      (IMAGINARY(fft_v2[i]) EQ 0.)) THEN BEGIN
        phase[i] = 0.
      ENDIF ELSE BEGIN
        phase[i] = ATAN(IMAGINARY(fft_v2[i]),FLOAT(fft_v2[i]))
      ENDELSE
    ENDELSE
  ENDELSE
ENDFOR

;define plot
LOADCT, 27
SET_PLOT, 'mac'
WINDOW, 0, XSIZE=650, YSIZE=500
;SET_PLOT, 'TEK'
;DEVICE, /tek4100
;!P.FONT = 0
;SET_PLOT, 'PS'
;LOADCT, 0
;DEVICE, FILENAME=plot_file_path+plot_filename, /PORTRAIT, $
;    XOFFSET=1.5, XSIZE=19.05, YOFFSET=1.35, YSIZE=12.7
;DEVICE, FILENAME=plot_file_path+plot_filename, /LANDSCAPE

; Defined the plot titles
bird = 589
;t_fun = 'b589_t_fun_p3_v1.dat'
;plot_title = 'Bird'+STRCOMPRESS(STRING(bird))+$
;    ' transfer function (b589_t_fun_p3ave.dat)'
x_title='frequency (mHz)'
y_title1= 'resp/rein'
y_title2= 'phase angle'
!mtitle = plot_title
!xtitle = x_title
!ytitle = y_title1
;!Y.MINOR = -1

; Make top panel plot
!MTITLE = plot_title
!XTITLE = '   '
!YTITLE = y_title1
!p.position = [0.15, 0.5, 0.95, 0.9]
PLOT, freqs_v2, ABS(fft_v2), XRANGE=[freqs_v2[0],freqs_v2[255]], XSTYLE=5
;PLOT, freqs_v2, ABS(fft_v2), XRANGE=[freqs_v2[0],freqs_v2[255]], XSTYLE=5, /YLOG
AXIS, XAXIS=0, XTICKNAME=REPLICATE(' ',30), XRANGE=[freqs_v2[0],freqs_v2[255]], XSTYLE = 1
AXIS, XAXIS=1, XTICKNAME=REPLICATE(' ',30), XRANGE=[freqs_v2[0],freqs_v2[255]], XSTYLE = 1

; Make bottom panel plot
plot_title = '  '
!MTITLE = plot_title
!XTITLE = x_title
!YTITLE = y_title2
!p.position = [0.15, 0.1, 0.95, 0.5]
PLOT, milli_freqs_v2, phase, XRANGE=[milli_freqs_v2[0],milli_freqs_v2[255]], XSTYLE=1, $
    YRANGE = [-!PI,!PI], YSTYLE=1, /NOERASE
AXIS, XAXIS=0, XTICKNAME=REPLICATE(' ',30), $
    XRANGE=[milli_freqs_v2[0],milli_freqs_v2[255]], XSTYLE = 1
!XTITLE = '   '
AXIS, XAXIS=1, XTICKNAME=REPLICATE(' ',30), $
    XRANGE=[milli_freqs_v2[0],milli_freqs_v2[255]], XSTYLE=1

;DEVICE, /CLOSE

;SET_PLOT, 'mac'

END