PRO PMTRAS_DBASE_MANAGER, SAVEFILE=savefile, TIME_RANGE=time_range, 				$
			DUMPMODE=dumpmode, PLOTMODE=plotmode, PHOTOMETRYMODE=photometrymode,  	$
			GAPMODE=gapmode, SPOTMODE=spotmode, DB_INTERPOLATE=db_interpolate,		$
			FITMODE=fitmode, PRED_TIMES=pred_times, PRED_ROLL=pred_roll,			$
			PLOT_FIT=plot_fit, PHOTO_RESID=photo_resid
;
; Multipurpose program to evaluate output of PMTRAS_DBASE_GENERATOR
;
; CHOICE OF DATABASE INPUT:
;	Default input is the most recent hsi_roll_db_master_yyyymmdd_hhmm.dat file.
; 	If SAVEFILE is specified, then the .SAV file (created by PMTRAS_DBASE_GENERATOR) is used instead.
;	If specified, time_range limits the range of times to be reviewed eg '2002-apr-21 '+['09:44','11:20']
;		(default = entire file)
;
; CHOICE OF ANALYSIS MODES:
; 	Default: Reviews file and notes abnormal data points by appending to PMTRAS_MGR_LOG.txt
; 	If DUMPMODE is set, only generates ascii display of PMTRAS database in structure, dbsolution
; 		with output written to PMTRAS_DBASE_DUMP.out
;	If PLOTMODE is set, only generates summary plots of PMTRAS database.
;	If PHOTOMETRYMODE is set, generates photometry-based output.
;	If GAPMODE is set, generates a file, ROLL_DBASE_GAPS.txt containing a compact list of gaps
;	IF SPOTMODE is set, generates a file, PMTRAS_MGR_SPOT.txt containing a list of values at evenly spaced intervals.
;	IF DB_INTERPOLATE is set to 1, 2 or 3 database is uses linear, quadratic or cubic interpolation
;			where feasible and generates a revised database.
;			Default = no interpolation. Using /DB_INTERPOLATION defaults to linear.
;	IF FITMODE is set, a multiparameter fit is used to fit roll angles.
;	PRED_TIMES = optional input array of sctimes at which the multiparameter fit is used to predict roll angles.
;	PRED_ROLL  = corresponding array of predicted roll angles (output).
;	PLOT_FIT   = 2-element vector containing the linear fit parameters, PA at time_range[0] and omega(radians/s).
;					A plot will be generated which shows the difference between roll data and the linear fit.
;
; Set bits in plotmode value select individual plots:
;	Bit 0 = 1 ==> AV_OFFSET vs time						(eg PLOTMODE = 1)
;	Bit 1 = 1 ==> AV_OFFSET histogram					(eg PLOTMODE = 2)
;	Bit 2 = 1 ==> Orbit-to-orbit phase shift vs time 	(eg PLOTMODE = 4)
;	Bit 3 = 1 ==> Roll_period vs time					(eg PLOTMODE = 8)
;	Bit 4 = 1 ==> Est_error vs time						(eg PLOTMODE = 16)
;	Bit 5 = 1 ==> Phase residual vs time				(eg PLOTMODE = 32)
;
; Modules used:
;	hsi_roll_db_full
;	hsi_any2sctime
;	PMTRAS_DBASE_DUMP
;	PMTRAS_DBASE_FIT
;	PMTRAS_DBASE_GAPS
;	PMTRAS_DBASE_PHOTOMETRY
;	PMTRAS_DBASE_PLOT
;	PMTRAS_DBASE_PREDICT
;	PMTRAS_DBASE_SPOTCHECK
;
; 14-nov-02 gh	Initial version checks for data gaps, no solution, and phase discontinuities.
; 15-nov-02 gh	Refine handling of data gaps.
;				Add provision to write output to log file.
; 11-Jan-03 gh	Add FILENAME keyword.
; 12-Jan-03 gh	Include sign in printed value of phase discontinuity
; 14-Jan-03 gh	Change default input to be most recent master data base file.
;				Change FILENAME keyword to SAVEFILE.
;				Change limit variable in main loop to LONG to support files > 48 days.
;  5-May-03 gh	Add TIME_RANGE keyword.
;				Add DUMPMODE keyword
; 10-Jul-03 gh	Add PLOTMODE keyword.
; 11-Jul-03 gh	Use plotmode as an argument to convey choice of plots.
; 16-Jul-03 gh	Debug inoperative SAVEFILE input mode.
; 23-JUl-03 gh	Add fair & uncertain lines to output log.
;				Set minimum quality to 'good' for discontinuity checks.
; 25-Jul-03 gh	Internalize call to ANYTIM for time_range argument.
;  1-Aug-03 gh	Add PHOTOMETRYMODE keyword.
;  9-Aug-03 gh	Add GAPMODE keyword.
; 11-Aug-03 gh	Add SPOTMODE keyword
; 13-Aug-03 gh	Add DB_INTERPOLATE keyword.
; 21-Aug-03 gh	Use successive linear/quadratic/cubic interpolation
;				Exploit /CREATE switch to write interpolated database
;  8-Mar-04 gh	Generate MESSAGE when roll database cannot be found.
; 15-Mar-O4 gh	Better implementation of 8-Mar change.
; 14-Apr-04 gh	Add FITMODE keyword.
; 15-Apr-04 gh	Add PRED_TIMES and PRED_ROLL keywords.
; 24-JUn-04 gh	Add photo_resid keyword to return average scatter in av_offset
;
; Define some constants and initialize counters.


pi64 = !PI * 64.													; used in extrapolating phases by 32s
PI63 = !PI * 63.
pi2	 = !PI * 2.
phaselim = 0.003													; phase matching tolerance (10 arcmin)
ngap 	= 0
ndisc	= 0
nosol	= 0
nfairsol = 0
nbadsol	= 0
;
; Open files and restrict input structure to requested time, if necessary.

IF N_ELEMENTS(savefile) NE 0 THEN RESTORE, savefile	ELSE BEGIN 		; restore dbsolution structure from .SAV file
	opmt = hsi_roll_db_full()
	dbsolution = opmt->getdata()									; restore dbsolution structure from master_db file
	IF N_ELEMENTS(dbsolution) EQ -1 THEN MESSAGE, 'PMTRAS_DBASE_MANAGER could not find the roll database.'
	obj_destroy, opmt
ENDELSE
IF N_ELEMENTS(time_range) NE 0 THEN BEGIN							; Just keep points within specified time_range.
	sctime_range 	= hsi_any2sctime(ANYTIM(time_range))
	inrange			= WHERE(dbsolution.sctime GE sctime_range[0].seconds AND 		$
							dbsolution.sctime LE sctime_range[1].seconds, npt)
	IF npt EQ 0 THEN STOP, 'No points within time_range.'
dbsolution		= dbsolution(inrange)
ENDIF
;
; In dump mode, use a separate procedure to generate output.
IF KEYWORD_SET(dumpmode) NE 0 THEN BEGIN
	PMTRAS_DBASE_DUMP, dbsolution
	RETURN
ENDIF
;
; In plot mode, use a separate procedure to generate output.
IF KEYWORD_SET(plotmode) NE 0 THEN BEGIN
	PMTRAS_DBASE_PLOT, dbsolution, plotmode
	RETURN
ENDIF
;
; In gap mode, use a separate procedure to generate output.
IF KEYWORD_SET(gapmode) NE 0 THEN BEGIN
	PMTRAS_DBASE_GAPS, dbsolution, gapmode
	RETURN
ENDIF
;
; In photometry mode, use a separate procedure to generate output.
IF KEYWORD_SET(photometrymode) NE 0 THEN BEGIN
	PMTRAS_DBASE_PHOTOMETRY, dbsolution, photometrymode, AVOFF_SCAT=photo_resid
	RETURN
ENDIF
;
; In spot mode, use a separate procedure to generate output.
IF KEYWORD_SET(spotmode) NE 0 THEN BEGIN
	PMTRAS_DBASE_SPOTCHECK, dbsolution
	RETURN
ENDIF
;
; Interpolate mode interpolates where feasible and generates a revised database if necessary.
IF KEYWORD_SET(db_interpolate) NE 0 THEN BEGIN
	PRINT, 'Gap statistics:         total number  short gaps  well-bracketed  filled'
	IF db_interpolate EQ 1 THEN BEGIN
		PMTRAS_DBASE_INTERPOLATE, dbsolution, 1, gapcounts, return_flag			; linear interpolation
		PRINT, gapcounts, FORMAT="('Linear interpolation:   ',  4I12)"
	ENDIF
	IF db_interpolate EQ 2 THEN BEGIN
		PMTRAS_DBASE_INTERPOLATE, dbsolution, 2, gapcounts, return_flag			; quadratic interpolation
		PRINT, gapcounts, FORMAT="('Quadratic interpolation:',  4I12)"
	ENDIF
	IF db_interpolate EQ 3 THEN BEGIN
		PMTRAS_DBASE_INTERPOLATE, dbsolution, 3, gapcounts, return_flag			; cubic interpolation
		PRINT, gapcounts, FORMAT="('Cubic interpolation:    ',  4I12)"
	ENDIF
	IF (return_flag NE 0) THEN BEGIN
	PRINT, 'WRITING NEW MASTER'
		opmt = hsi_roll_db_full()
		opmt->update, dbsolution=dbsolution, /CREATE			; create new master_db file
		obj_destroy, opmt
	ENDIF
	RETURN
ENDIF
;
; In fit mode, use a separate procedure to generate output.
IF KEYWORD_SET(fitmode) NE 0 THEN BEGIN
	PMTRAS_DBASE_FIT, dbsolution
	RETURN
ENDIF
;
; In prediction mode, fitmode is invoked first, then a separate function is called to generate output.
IF N_ELEMENTS(pred_times) NE 0 THEN BEGIN
	PMTRAS_DBASE_FIT, dbsolution, fitparm
	pred_roll = PMTRAS_DBASE_PREDICT(pred_times-dbsolution[0].sctime, fitparm)	; note that times are now relative!
	pred_roll = pred_roll* !DTOR				; convert from degrees to radians
	RETURN
ENDIF
;
; In plot_fit mode, generate plot to compare data to a linear fit to pred_roll.
IF N_ELEMENTS(plot_fit) EQ 2 THEN BEGIN
	omega		= plot_fit[1]												; radians/second
	pa_0		= plot_fit[0]												; radians
	t0			= DOUBLE(sctime_range[0].seconds)							; seconds
	temp		= pb0r(anytim(time_range[0]))								; p,b0,r in deg, deg, arcmin
	solar_pa	= temp[0] * !DTOR											; in radians
	ok 			= WHERE(dbsolution.roll_quality GE 192)
	trelok		= dbsolution[ok].sctime - t0
	rollok		= dbsolution[ok].roll_phase/1.D6 - solar_pa				; converting from microradians to radians
	fit_roll 	= pa_0 + omega*trelok									; linear fit
	residual	= (rollok - fit_roll) / (2*!DPI) +0.5				; (observed angle - linear fit) in rotations
	resid		= (residual - FLOOR(residual) -0.5) * 360.							; in range [-180, +180]
	t0_ref		= hsi_sctime2any({seconds: dbsolution[0].sctime, bmicro: 0L})
	UTPLOT, trelok, resid, t0_ref, YTITLE='ACTUAL ROLL ANGLE - LINEAR FIT (degrees)', PSYM=6
	OPLOT, trelok, fit_roll-fit_roll
	RETURN
ENDIF
;
; Create time display strings.
now = SYSTIME()
current_time = STRMID(now,8,2)+'-'+STRMID(now,4,3)+'-'+STRMID(now,20,4)+STRMID(now,10,6) ; 'yyyy-mmm-dd hh:mm'
OPENU, dbmlogfile, 'PMTRAS_MGR_LOG.txt', /GET_LUN, /APPEND
npt 					= N_ELEMENTS(dbsolution)
sc_event_time 			= {seconds: 0L, bmicro: 0L}						; define a proper s/c time structure
sc_event_time 			= REPLICATE(sc_event_time, npt)
sc_event_time.seconds 	= dbsolution.sctime
showtime				= hsi_sctime2any(sc_event_time, /STIME)			; = array of time-display strings.
PRINTF, dbmlogfile
PRINTF, dbmlogfile, current_time, showtime[0], FORMAT="(A17, 3X, A20, 8X, 'START')"
;
; Begin loop over data points in dbsolution (excluding the last one).
dphtab = FLTARR(npt-1)
MODTIME = fltarr(npt-1)
FOR n=0L, npt-2 DO BEGIN										; Begin looping over data points.
	IF dbsolution[n].roll_quality LT 64 THEN BEGIN											; if no solution
		PRINTF, dbmlogfile, current_time, showtime[n], FORMAT="(A17, 3X, A20, 8X, 'No solution.')"
		nosol = nosol + 1
	ENDIF ELSE BEGIN
		IF dbsolution[n].roll_quality LT 128 THEN BEGIN
			PRINTF, dbmlogfile, current_time, showtime[n], FORMAT="(A17, 3X, A20, 8X, 'Bad solution.')"
			nbadsol = nbadsol + 1
		ENDIF ELSE BEGIN
			IF dbsolution[n].roll_quality LT 192 THEN BEGIN
				PRINTF, dbmlogfile, current_time, showtime[n], FORMAT="(A17, 3X, A20, 8X, 'Fair solution.')"
				nfairsol = nfairsol + 1
			ENDIF ELSE BEGIN
				dt = dbsolution[n+1].sctime - dbsolution[n].sctime
				IF dt GT 64.1 OR dt LE 63.9 THEN BEGIN																; if gap
					PRINTF, dbmlogfile, current_time, showtime[n], dt, FORMAT = "(A17, 3X, A20, F7.0, ' second gap.')"
					ngap = ngap + 1
				ENDIF ELSE BEGIN
					IF dbsolution[n+1].roll_quality LT 192 THEN CONTINUE			; skip if next solution not good
					phfwd = dbsolution[n].roll_phase  *1.E-6 + pi64 / dbsolution[n].roll_period
					phbak = dbsolution[n+1].roll_phase*1.E-6 - pi64 / dbsolution[n+1].roll_period
					dph	= (phfwd - phbak + pi63) MOD pi2 - !PI
					IF ABS(dph) GT phaselim THEN BEGIN													; if discontinuity
						PRINTF, dbmlogfile, current_time, showtime[n], dph, 	$
								FORMAT = '(A17, 3X, A20, F7.3, " radian phase discontinuity.")'
						ndisc = ndisc + 1
					ENDIF
				ENDELSE
			ENDELSE
		ENDELSE
	ENDELSE
ENDFOR
;
; Special handling of last data point
IF dbsolution[npt-1].roll_quality LT 64 THEN BEGIN											; if no solution
	PRINTF, dbmlogfile,  current_time, showtime[npt-1], FORMAT="(A17, 3X, A20, 8X, 'No solution.')"
	nosol = nosol + 1
ENDIF ELSE BEGIN
	IF dbsolution[npt-1].roll_quality LT 128 THEN BEGIN											; if bad solution
		PRINTF, dbmlogfile,  current_time, showtime[npt-1], FORMAT="(A17, 3X, A20, 8X, 'Bad solution.')"
		nobadsol = nobadsol + 1
	ENDIF ELSE BEGIN
		IF dbsolution[npt-1].roll_quality LT 64 THEN BEGIN											; if fair solution
			PRINTF, dbmlogfile,  current_time, showtime[npt-1], FORMAT="(A17, 3X, A20, 8X, 'Fair solution.')"
			nfairsol = nfairsol + 1
		ENDIF
	ENDELSE
ENDELSE
PRINTF, dbmlogfile,  current_time, showtime[npt-1], FORMAT="(A17,3X, A20, 8X, 'END')"
PRINTF, dbmlogfile
FREE_LUN, dbmlogfile
PRINT, 'PMTRAS_DBASE_MANAGER done.'
END