;+
; Project     : ovsa
;
; Name        : OVSA__DEFINE
;
; Purpose     : Define an ovsa object
;
; Category    : OVSA analysis
;
; Explanation : READ,  PLOT and HELP methods for OVSA FITS
;
; Syntax      : This procedure is invoked when a new ovsa object is
;               created
;
; Examples    : Plot a total power light curve at 5.8 GHz from antenna 1 
;               IDL> a = OBJ_NEW( 'ovsa' )
;               IDL> a->READ, 'ovsa.fits', data
;               IDL> a->PLOT, 1, 5.8, /lightcurve
;
; Inputs      : 'ovsa' = object classname
;              
; Outputs     : Object with methods: READ,  PLOT and HELP
;
; Keywords    : None
;                   
; History     : Written 8-Dec-1999 ptg
;               Added movie keyword 4-April-2000 ptg 
;               Added fit keyword 22-May-2000 ptg
;               Added MESH keyword and checked code 21-Nov-2000 ptg
;
; Contact     : ptg@bbso.njit.edu (Peter Gallagher,  NJIT)
;-

;----------------------------------------------------
; Define help method for ovsa class

pro ovsa::help
 
  iprint, obj_methods( 'ovsa' )

return
end

;----------------------------------------------------
; Define read method for ovsa class

pro ovsa::read, filename, data, quiet = quiet


; Read ovsa FITS file
 
  read_ovsa_fits, filename, data
  
  
; Print information about observations.
  if ( keyword_set( quiet ) eq 0 ) then begin  
    print,  ''
    print,  '*********** Observation Information Follows *********'
    print,  ''
    print,  data.header( 12 )
    print,  data.header( 14 )
    print,  data.header( 15 )
    print,  data.header( 16 )
    print,  data.header( 17 )
    print,  data.header( 18 )
    print,  data.header( 19 )
    print,  data.header( 20 )
    print,  'ANTENNA =             ' + data.feedinfo( 0 ) + $
            '/Antenna number and polarization'
    for i = 1,  n_elements( data.feedinfo ) - 1 do $
    print,  '                      ' + data.feedinfo( i ) + '/'
    print,  ''
    print,  '********************** End ***************************'
    print, ' '
  endif

; Associate pointers with heap arrays
  
  self.header = ptr_new( data.header )
  self.avg    = ptr_new( data.avg )
  self.tavg   = ptr_new( data.tavg )
  self.f      = ptr_new( data.f )

return
end


;----------------------------------------------------
; Define plot method for ovsa class
; NOTE: The argument gen can be the frequency or the time required.
; 

pro ovsa::plot, ant, gen, tmin, tmax, dynamic = dynamic, $
                lightcurve = lightcurve, spectrum = spectrum, $
	        interactive = interactive, save = save, $
		movie = movie, fit = fit, mesh = mesh, $
		lin = lin, lcp = lcp, rcp = rcp, $
		background = background, yrange = yrange, $
                quiet = quiet, title = title


; Check if the antenna number and polarization choices are correct
  
  case ant of
    1:begin 
        if ( keyword_set( lin ) ) then begin   ; Ant. 1 linear
	  channel = 0 & pol = 'LIN'
	endif 
        if ( keyword_set( lcp ) ) then begin   ; Ant. 1 left circular
	  channel = 1 & pol = 'LCP'  
	endif
      end
    2:begin 
        if ( keyword_set( rcp ) ) then begin   ; Ant. 2 right circular
	  channel = 2 & pol = 'RCP'  
	endif  
        if ( keyword_set( lin ) ) then begin   ; Ant. 2 linear
	  channel = 3 & pol = 'LIN'  
	endif  
      end
    3:                                         ; Ant. 3 not available
    4:begin
        if ( keyword_set( lin ) ) then begin   ; Ant. 4 LIN
          channel = 4 & pol = 'LIN'  
	endif 
      end	 
    5:begin 
        if ( keyword_set( lin ) ) then begin   ; Ant. 5 LIN
	  channel = 5 & pol = 'LIN'  
	endif
      end	
    6:begin  
        if ( keyword_set( lin ) ) then begin   ; Ant. 6 LIN
	  channel = 6 & pol = 'LIN'  
	endif
      end	  
    7:begin
        if ( keyword_set( lin ) ) then begin   ; Ant. 7 LIN 
	  channel = 7 & pol = 'LIN'  
	endif
      end	  
  endcase
  
  if ( n_elements( channel ) eq 0 ) then begin
    print, ''
    print, 'OVSA> Invalid Antenna Number or Polarization Choice'
    print, ''
    print, 'OVSA> Use:  Antenna 1 with keyword lin or lcp'
    print, '            Antenna 2 with keyword lin or rcp'
    print, '            Antenna 4 with keyword lin'
    print, '            Antenna 5 with keyword lin'
    print, '            Antenna 6 with keyword lin'
    print, '            Antenna 7 with keyword lin'
    goto, jump
  endif


; Extract frequency and time array

  t = (*self.tavg) / 1000.         ; Convert time to seconds
  f = *self.f    
  t_start = strmid( (*self.header)( 14 ), 11, 23 )
  t_end   = strmid( (*self.header)( 15 ), 11, 23 )


; If times not input as arguments, then define start
; and end times
  
  if ( n_params() le 2 ) then begin
      	tmin = strmid( t_start, 11, 12 ) 
        tmax = strmid( t_end,   11, 12 )
  endif

        
; Regrid data if the mesh keyword is set and do background subtraction
; if the background keyword is set.

  iold = reform( (*self.avg)( channel,  *,  * ) )    	
  
  if ( keyword_set( mesh ) ) then begin  
    ovsa_mesh,  iold,  f,  t,  inew,  fnew,  tnew
    iold = reverse( rotate( iold,  3 ), 2 )
    inew = reverse( rotate( inew,  3 ), 2 )
    fnew = exp( fnew )
    if ( keyword_set( background ) ) then inew = ovsa_bg( inew )
  endif else begin
    iold( where ( finite( iold ) eq 0 ) ) = 0.001  ; Set NaNs to 0.001!
    iold = reverse( rotate( iold, 3 ), 2 )         ; ... NEEDS TO BE CHANGED!
    inew = iold
    if ( keyword_set( background ) ) then inew = ovsa_bg( inew )
    fnew = f
    tnew = t
  endelse  
  
  
; Run OVSA_INTERACT if the interactive mode is selected
  
  if ( keyword_set( spectrum ) and keyword_set( interactive ) ) then begin
     ovsa_interact, f, iold, inew, tnew, fnew, ant, t_start, /spectrum
     goto, jump
  endif  

  if ( keyword_set( lightcurve ) and keyword_set( interactive ) ) then begin
     ovsa_interact, f, iold, inew, tnew, fnew, ant, t_start, /lightcurve
     goto, jump
  endif
  
  if ( (n_params() gt 1) and keyword_set( dynamic ) ) then begin 
    tmax = tmin 
    tmin = gen
  endif 


; Change time format from hh:mm to hh:mm:ss if needed.

  if ( strlen( tmin ) le 5 ) then tmin = tmin + ':00'
  if ( strlen( tmax ) le 5 ) then tmax = tmax + ':00'


; Locate pixels in dynamic spectrum corresponding to the requested 
; minimum and maximum times.

  tsmin = hms2sec( strmid( tmin, 0, 2 ) + $
                   strmid( tmin, 3, 2 ) + $
	           strmid( tmin, 6, 4 ) ) - tnew( 0 )
		   
  tsmax = hms2sec( strmid( tmax, 0, 2 ) + $
                   strmid( tmax, 3, 2 ) + $
		   strmid( tmax, 6, 4 )) - tnew(0)
		   
  t = tnew - tnew(0)
  
  t_low  = where( ( t - tsmin )^2 eq min( ( t - tsmin )^2 ) )
  t_high = where( ( t - tsmax )^2 eq min( ( t - tsmax )^2 ) )
  
  if ( (fix( tsmin ) lt fix( t( 0 ) ) ) or $
       (fix( tsmax ) gt fix( t( n_elements( t ) - 1 ) ) ) ) then begin
    print, ''
    print, 'OVSA> Times must be between ' + $
                  strmid( t_start, 11, 8) + $
		  ' and ' + $
		  strmid( t_end, 11, 8)
    print, ''
  endif
  
  
; Plot lightcurve.
          
  if keyword_set( lightcurve ) then begin
   
    freq = gen

    ; Find frequency closest to requested freq.	 
       
    freqn_index = where( ( f - freq )^2 eq min( ( f - freq )^2 ) )
    freqn_obs = f( freqn_index )
    freqn_index = freqn_index( 0 )
    freqn_obs = strmid( arr2str( freqn_obs( 0 ), /trim ), 0, 4 )  
  			 ; If the requested frequency lies directly
			 ; between 2 observed frequencies,  choose the
     			 ; lower.
    
    if ( keyword_set( quiet ) eq 0 ) then begin
      print, ''
      print, 'OVSA> Requested Frequency:	     ' + $
              arr2str( freq, /tr) + ' GHz'
      print, 'OVSA> Closest Observed Frequency: ' + $
              arr2str( freqn_obs, /tr ) + ' GHz'
      print, ''
    endif

    if ( keyword_set( yrange ) eq 0 ) then $
         yrange = [ 0, max( reform( inew( *, freqn_index ) ) ) ]
    if ( keyword_set( title ) eq 0 ) then $
         title = 'Total Power For Antenna ' + arr2str( ant, /trim ) + $
                 ' (' + pol + ') at ' + freqn_obs + ' GHz' 

    utplot, tnew - min( tnew ), reform( inew( *, freqn_index ) ), t_start, $
	    title = title, /xs, $ 
	    ytitle = 'Flux/SFU', psym = 10, $
	    timerange = [strmid( t_start, 0, 11 ) + tmin, $
	                 strmid( t_start, 0, 11 ) + tmax ], $
	    yrange = yrange
		 
    print, 'OVSA> Time range available: ' + strmid( t_start, 11, 8 ) + $
                  ' to ' + strmid( t_end, 11, 8 )
    print, 'OVSA> Time range plotted:   ' + strmid( tmin, 0, 8 ) + $
                  ' to ' + strmid( tmax, 0, 8 )
    print, ''
	 
    if keyword_set( save ) then begin
      ovsa_data = { t_sec  : tnew - min( tnew ), $
                    power  : reform( inew( *, freqn_index ) ), $
		    t_base : t_start, $
     		    t_start: strmid( t_start, 0, 11 ) + tmin, $
	  	    t_end  : strmid( t_start, 0, 11 ) + tmax, $
      		    freq   : arr2str( freqn_obs, /tr ) }
		   
      save, file='ovsa.sav', ovsa_data
      print, 'OVSA> Data written to ovsa.sav'
      print, ''
    endif

  endif
  

; Plot dynamic spectrum.
 
  if keyword_set( dynamic ) then begin

    if keyword_set( movie ) then begin
      ovsa_movie, tmin, tmax, t_start, t_end, t_low, t_high, $
   		  fnew, inew, tnew, ant, iold, f, t
       goto, jump
    endif
   
    pm = !p.multi( 0 )

    if ( keyword_set( title ) eq 0 ) then $
         title = 'Dynamic Spectrum For Antenna ' + arr2str( ant, /trim ) + $
                 ' (' + pol + ')'

    plot_image, reverse( inew( t_low( 0 ) : t_high( 0 ), * ), 2 ), $
   	        /nosquare, xticks = 1, yticks = 1, $
   	        xtickname = [" ", " "], ytickname = [" ", " "], $
   	        title = title, ytitle = 'Frequency/GHz'

    am = !p.multi( 0 )
    !p.multi( 0 ) = pm
    
    utplot, tnew - min( tnew ), reverse( fnew ), t_start, /noerase, $
    	    /nodata, /xs, /ys, yr=[ max( fnew ), min( fnew ) ], $
    	    timerange=[ strmid( t_start, 0, 11 ) + tmin, $
   	        	strmid( t_start, 0, 11 ) + tmax ]
    !p.multi( 0 ) = am
    
    if ( keyword_set( quiet ) eq 0 ) then begin	
      print, ''
      print, 'OVSA> Time range available: ' + $
                    strmid( t_start, 11, 8 ) + ' and ' + $
   		  strmid(   t_end, 11, 8 )
      print, 'OVSA> Time range plotted:   ' + $
                    strmid( tmin, 0, 8 ) + ' and ' + $
                    strmid( tmax, 0, 8 )
      print, ''
    endif

    if keyword_set( save ) then begin
      ovsa_data = { image_new : reverse( inew( t_low( 0 ) : t_high( 0 ), * ), 2 ), $
   		    t_base    : t_start, $
   		    t_sec     : tnew - min( tnew ), $
   		    f_ghz     : fnew, $
   		    image_old : iold, $
   		    f_old     : f, $
   		    t_start   : strmid( t_start, 0, 11 ) + tmin, $
   		    t_end     : strmid( t_start, 0, 11 ) + tmax }
      save, file='ovsa.sav', ovsa_data
      print, 'OVSA> Data written to ovsa.sav'
      print, ''
    endif
   		   
  endif


; Plot spectrum at selected time

  if keyword_set( spectrum ) then begin
    
    time = gen
    if ( strlen( time ) le 5 ) then time = time + ':00'
      
     ;Locate pixels in dynamic spectrum corresponding to the requested 
     ;minimum and maximum times
      
    timesec = hms2sec( strmid( time, 0, 2 ) + $
    	               strmid( time, 3, 2 ) + $
                       strmid( time, 6, 2 ) ) - tnew( 0 )
    tsec = where( ( t - timesec )^2 eq min( ( t - timesec )^2 ) )
    tsec = tsec( 0 )
    
    if ( ( fix( timesec ) lt fix( t( 0 ) ) ) or $
         ( fix( timesec ) gt fix( t( n_elements( t ) - 1 ) ) ) ) then begin
      print, ''
      print, 'OVSA> Time must be between ' + $
                    strmid( t_start, 11, 8 ) + ' and '+ $
		    strmid(   t_end, 11, 8 )
      print, ''
      timesec = t( 0 )
    endif else begin
      print, ''
      print, 'OVSA> Time range available: ' + $
                    strmid( t_start, 11, 8 ) + ' and ' + $
		    strmid( t_end,   11, 8 )
      print, ''
      circle_sym, /fill
      
      if ( keyword_set( yrange ) eq 0 ) then $
           yrange = [ min( iold( tsec, * ) ), max( iold( tsec, * ) ) ]
      if ( keyword_set( title ) eq 0 ) then $
           title = 'Spectrum For Antenna ' + arr2str( ant, /trim ) + $
                        ' (' + pol + ') at ' + arr2str( time, /trim ) + ' UT'
 
      plot_oo, f, inew( tsec, * ), psym = 8, /xs, $
    	       xtitle = 'Frequency/GHz', $
               ytitle = 'Total Power/SFU', $
               title  = title, $
	       yrange = yrange, xrange = [ 1, 20]
      oplot, f, inew( tsec, * ), psym = -1, symsize = 2 

      if keyword_set( fit ) then begin
	yfit = ovsa_fit( f, inew( tsec, * ), /print, /inter, yrange = yrange )
        oplot, yfit( 0, * ), yfit( 1, *), thick = 2
      endif

    endelse
	             	    
    if keyword_set( save ) then begin
      ovsa_data = { f_ghz : f, $
	            power : reform( iold( tsec, * ) ), $
		    time  : arr2str( time, /trim ), $
		    yfit   : yfit }
      save, file='ovsa.sav', ovsa_data
      print, ''
      print, 'OVSA> Data written to ovsa.sav'
      print, ''
    endif
 	              
  endif
      
  jump:
  	 
return
end


; Set up ovsa object class

pro ovsa__define
   	      
  ovsa_struct = { ovsa, header : ptr_new(), $
                        avg    : ptr_new(), $
                        tavg   : ptr_new(), $
		        f      : ptr_new() }
return
end