;+
; NAME:
;	xhcs
; PURPOSE:
;	 Read a file that lists a set of Nick's WSO or NSO files.
;	 Run X's code on them, and create source surface maps
; CALLING SEQUENCE:
	PRO xhcs, obs, $
		destination		= destination	, $
		list			= list			, $
		source			= source		, $
		remote_prefix	= remote_prefix	, $
		remote_postfix	= remote_postfix, $
		local_prefix	= local_prefix	, $
		local_postfix	= local_postfix	, $
		gzip			= gzip			, $
		rss				= rss			, $
		rcp				= rcp			, $
		silent			= silent
; INPUTS:
;	obs			scalar; type: character; default: 'wso'
;					indentifies origin of magnetic field data
;					Currently supported:
;					wso_noaa	New Fits files from NOAA
;					wso_arge	Old .dat files from NOAA
;					nso_noaa	Old Fits files from NOAA
;					nso_ktpk	New Fits files from Kit Peak
; OPTIONAL INPUTS:
;	destination=destination
;				scalar; type: string; default: $TUB
;					directory where to write resulting source
;					surface magnetic field (the file name is the
;					same as for the file with the input photospheric
;					magnetic field.
;	list=list	scalar; type: string; default: none
;					name of file containing fully-qualified file names
;					for all photospheric synoptic maps to be processed.
;	source=source
;				scalar; type: string; default: none
;					used only if keyword list is not specified.
;					Name of single photospheric magnetic map to be
;					processed. 
; INCLUDE:
	@compile_opt.pro		; On error, return to caller
; CALLS:
;	InitVar, txt_read, SetFileSpec, GetFileSpec, brss_csss_scbf_
;	boost, wso_write, destroyvar, IsType, hide_env
; RESTRICTIONS:
;	The input should not contain any holes. This sometimes happens with
;	the WSO data received from NOAA when there are data gaps.
;	Xuepu's code seems to return an array of zeros when this happens.
; MODIFICATION HISTORY:
;	Jun 2002 (tjd, tdunn@ucsd.edu)
;	Spt 2002 (tjd, tdunn@ucsd.edu) removed the intermediary 
;		printing of files to /xwso/sync, which formerly 
;		worked as a checking procedure, and is now commented out.  
;		Made code flexible to changes in $dat dir.
;	NOV-2002, Paul Hick (UCSD/CASS)
;		Modified XuePus readwsosc so that it can directly read
;		Nicks files (it now calls wso_read).
;		Writing of source surface magnetic field is now done using wso_write.
;	OCT-2003, Paul Hick (UCSD/CASS)
;		Modified to intercept photospheric maps with data gaps. These are
;		now skipped. The modification required adding the keyword 'badfile'
;		to some of Xuepus original code (gh_wbrsc, readwsosc).
;	OCT-2005, Paul Hick (UCSD/CASS; pphick@ucsd.edu)
;		Added keyword remote_prefix to allow for a change in prefix from
;		input file to output file (used for NSO/SOLIS maps).
;-

InitVar, obs		, 'nso_ktpk'			; Observatory
InitVar, remote_prefix	, obs
InitVar, remote_postfix	, ''
InitVar, local_prefix	, obs
InitVar, local_postfix	, ''
InitVar, destination	, getenv('TUB')
InitVar, gzip		, /key
InitVar, rss		, 15.0		; Source surface distance (solar radii)
InitVar, rcp		,  1.6		; ?? Cusp distance (solar radii)
InitVar, silent		, 0

ghf = filepath(root=getenv('TUB'),'ghf.txt') ; Spherical harmonic expansion coeff

; Read master file

CASE IsType(list, /defined) OF
0: BEGIN							; No list file, use keyword source
	nfile = n_elements(source)
	IF nfile GT 0 THEN infile = source
END
1: BEGIN
	IF NOT txt_read(list, infile, nrec=nfile, /silent) THEN BEGIN
		message, /info, 'error reading master file '+hide_env(list)
		RETURN
	ENDIF
	IF silent LE 0 THEN message, /info, 'master file: '+hide_env(list)+' ('+strcompress(nfile,/rem)+' files)'
END
ENDCASE

IF IsType(infile, /undefined) THEN BEGIN
	message, /info, 'no photospheric maps specified'
	RETURN
ENDIF

indir = GetFileSpec(infile,upto='directory')
indir = indir[ uniq(indir,sort(indir)) ]

IF n_elements(indir) GT 1 THEN BEGIN
	message, /info, 'all photospheric maps must be in same directory'
	RETURN
ENDIF

indir = indir[0]
destination = GetFileSpec( filepath(root=destination, 'dummy.file'), upto='directory')

IF destination EQ indir THEN BEGIN
	message, /info, 'input and output directory must be different'
	REURN
ENDIF

name = GetFileSpec(infile, part='name',/strict)
IF local_prefix  NE remote_prefix  THEN name = local_prefix+strmid(name,strlen(remote_prefix))
IF local_postfix NE remote_postfix THEN name = strmid(name,0,strlen(name[0])-strlen(remote_postfix))+local_postfix
outfile = filepath(root=destination, name)+'.fts'

; Read the list of files flagged as bad by previous runs

badlist = filepath(root=getenv('DAT'),subdir='map','bad_'+obs+'_obs.txt')
IF (file_search(badlist))[0] NE '' THEN sts = txt_read(badlist, old_badfiles, /silent)
destroyvar, new_badfiles

FOR i=0,nfile-1 DO BEGIN

	destroyvar, header			; Safety belt, probably not necessary

	; Xuepus original readwsosc procedure reads only the first 72 longitudes from
	; the input file. The last record (same longitude as first, one rotation later)
	; is discarded.
	; This was changed so that now the first (oldest) longitude is discarded,
	; and the last (most recent) longitude is retained.
	;
	; cvar contains Carrington variables for all 72 longitudes.

	brss_csss_scbf_, infile[i], obs, ghf, rcp=rcp, rss=rss, brad=brad, cvar=cvar, header=header, badfile=badfile

	CASE badfile OF

	0: BEGIN

		IF silent LE 1 THEN message, /info, hide_env(infile[i])+' > '+hide_env(outfile[i])

		; Assume that brad is a 72x30 array with longitude (1st dim) and
		; sine(latitude) (2nd dim) increasing with array index in equal steps.

		boost, header, [	$
			'RSS='  +strcompress(rss,/rem)+' solar radii',	$
			'RCUSP='+strcompress(rcp,/rem)+' solar radii',	$
			'BMAP='+GetFileSpec(infile[i],from='name')   ,	$
			'DESCRIPT=Radial magnetic field at source surface from brss_csss_scbf',	$
			''	]

		destroyvar, obstime
		wso_write, outfile[i], brad, cvar=cvar, header=header, /fill, badfillvalue=0.0, /addfirst, obstime=obstime
		IF gzip THEN status = gzip_file(outfile[i],/force)

		IF IsType(obstime,/defined) THEN BEGIN
			openw, /get_lun, iu, filepath(root=getenv('DAT'),subdir='map','last_'+obs+'_obs.txt')
			printf, iu, obstime, '  ', GetFileSpec(out_file,from='name')
			free_lun, iu
		ENDIF
	END

	1: BEGIN

		badfile = GetFileSpec(infile[i],from='name')
		IF n_elements(old_badfiles) EQ 0 THEN BEGIN
			boost, new_badfiles, badfile
		ENDIF ELSE IF (where(old_badfiles EQ badfile))[0] EQ -1 THEN BEGIN
			boost, new_badfiles, badfile
		ENDIF

	END

	ENDCASE

ENDFOR

; Overwrite the badfiles list if new bad files were encountered in this run

n = n_elements(new_badfiles)
IF n NE 0 THEN BEGIN
	FOR i=0,n-1 DO message, /info, 'file skipped; probably bad data: '+new_badfiles[i]
	message, /info, '(file(s) added to '+hide_env(badlist)+')'
	boost, old_badfiles, new_badfiles
	openw, /get_lun, iu, badlist
	FOR i=0,n_elements(old_badfiles)-1 DO printf, iu, old_badfiles[i]
	free_lun, iu
ENDIF

RETURN  &  END