C+
C NAME:
C	ReadVIPSn8_n
C	ReadVIPSLOSCheckn_n
C PURPOSE:
C	Read IPS velocity data into arrays
C CATEGORY:
C	I/O
C CALLING SEQUENCE:
C	call ReadVIPSn8_n(iSys,iReadGen8_n,iProcessGen8_n,cWildVIPS,
C     &		nCar,JDCar,bAuto,XCtest1,XCtest2,NCoff,MJDref,
C     &		Radius,
C     &		XElow,XEhigh,iEorW,XRlim,
C     &		NLmax,NL,cbvfiles,SRCVG,SRC,SRCVGsav,SRCsav,
C     &		IYRF,IREC,IYRS,DOYS8,XDS,XLS,XLL,XDL,XCE,XE,XC,YL,VV,DVV,GG,DGG,FRE,SSS,
C     &		NSmax,IYRFsav,IRECsav,IYRSsav,DOYSsav8,XDSsav,XLSsav,XLLsav,
C     &		XDLsav,XCEsav,XEsav,XCsav,YLsav,VVsav,DVVsav,GGsav,DGGsav,FRsav,SSSsav)
C
C	call ReadVIPSLOSCheckn_n(XCbeg,XCend,NLOS,NINSIDE,DLOS)
C INPUTS:
C   ReadVIPSn8_n:
C	iSys		integer		integer indicating which IPS system is used (modified 8/30/2013 BVJ)
C					1 = UCSD
C					2 = Cambridge
C					3 = Ooty
C					4 = Nagoya
C					5 = KSWC
C					6 = MEXART
C					7 = Pushchino
C					8 = LOFAR
C					if all are negative, this is a g-level read
C	iReadGen8_n	external integer function
C	iProcessGen8_n	external integer function
C					external functions control data file access
C	cWildVIPS	character*(*)	wildcard used to locate data files;
C					the wildcard must contain the substring 
C					'%%'
C	nCar		integer		dimension of JDCar
C	JDCar(nCar)	integer		start times (Julian days) of Carrington rotations
C	bAuto		logical
C	XCtest1		real		modified Carrington variable for start of search
C	XCtest2		real		modified Carrington variable for end of search
C	NCoff		integer		JDCar(I) is the start of rotation NCoff+I
C	MJDref		real*8		if not equal 0, only data prior to
C					modified Julian day MJDref are used
C	Radius		real		reference distance (AU)
C
C	XElow		real
C	XEhigh		real
C	iEorW		integer
C	XRlim		real
C
C	NLmax		integer		max. # points read into arrays 
C					XC,YL,VV,DVV,GG,DGG XCsav,YLsav,VVsav,DVVsav,GGsav,DGGsav
C	NSmax		integer
C
C    ReadVIPSLOSCheckn_n:
C	XCbeg		real		modified Carrington variable for start of map
C	XCend		real		modified Carrington variable for end of map
C	NLOS		integer
C	NINSIDE		integer
C	DLOS
C OUTPUTS:
C	NL		integer		# valid data points in arrays XC,YL,VV,DVV,GG,DGG
C	SRC(NLmax)	character*8	Character Source identifier
C	IYRF(NLmax)	integer		year of file from which observation was extracted (should be same as IYRS)
C	IREC(NLmax)	integer		record number on file IYRF
C	DOYS8(NLmax)	real*8		time of observation: day of year (incl. fraction for time of day)
C	XDS (NLmax)	real		Sun-Earth distance
C	XLS (NLmax)	real		geocentric ecliptic longitude Sun
C	XLL (NLmax)	real		geocentric ecliptic lng(P)-lng(Sun) (deg)
C	XDL (NLmax)	real		geocentric ecliptic lat(P) (deg)
C	XCE (NLmax)	real		modified Carrington variable of sub-Earth point on Sun
C	XE  (NLmax)	real		elongation (deg) (>0: East of Sun; <0: West of Sun)
C	XC  (NLmax)	real		modified Carrington variable point P after traceback to
C					  heliocentric distance Radius at speed VV
C	YL(NLmax)	real		heliographic latitude point P
C	VV(NLmax)	real		IPS velocity (km/s)
C	DVV(NLmax)	real		IPS velocity error (km/s)
C	GG(NLmax)	real		IPS disturbance factor
C	DGG(NLmax)	real		IPS disturmance factor error
C	FRE(NLmax)	real		IPS radio observation frequency
C	SSS(NLmax)	real		IPS source size
C	NSmax		integer		# data points saved in XCsav,YLsav,VVsav
C	SRCsav(NSmax)	character	scratch arrays (used internally only),
C	IYRFsav(NSmax)	integer		
C	IRECsav(NSmax)	integer
C	IYRSsav(NSmax)	integer
C	DOYSsav(NSmax)	real*8
C	XDSsav (NSmax)	real
C	XLLsav (NSmax)	real
C	XDLsav (NSmax)	real
C	XCEsav (NSmax)	real
C	XEsav  (NSmax)	real
C	XCsav  (NSmax)	real
C	YLsav  (NSmax)	real
C	VVsav  (NSmax)	real
C	DVVsav (NSmax)	real
C	GGsav  (NSmax)	real
C	DGGsav (NSmax)	real
C	FRsav  (NSmax)	real
C	SSSsav (NSmax)	real
C RESTRICTIONS:
C	Data files are read using a logical unit number assigned by iGetLun.
C	Since a data file may be left open on return, there is a potential
C	danger if the same unit is used somewhere. Assign unit numbers
C	with iGetLun to make sure this doesn't happen.
C FUNCTIONS/SUBROUTINES:
C	bOpenFile, itrim, SunNewcomb, iSetFileSpec, iGetFileSpec, Julian
C EXTERNAL:
C INCLUDE:
C	include		'filparts.h'
C	include		'sun.h'
C COMMON BLOCKS:
C SIDE EFFECTS:
C PROCEDURE:
C >	All points in the range [XCtest1,XCtest2] are extracted from the data files
C >!!!!	The modified Carrington variables are in units offset by NCoff from
C	the regular Carrington rotation number (i.e. in units which can be
C	used as index to JDCar)
C >	Two external user-defined are needed to read the data files.
C	External functions exist for reading the Nagoya and the UCSD velocity
C	IPS data:  iReadGen8_n   , iProcessGen8_n
C
C	These are programmed in the form:
C
C	function iReadNagoyan8_n(iSys,iU,iFirst)
C	... ( read a record from the data file ) ...
C	return
C
C	Input values to iReadGen8_n
C	iU	integer		logical unit number assigned by bOpenFile
C	iFirst	integer		1 : initializes record counter
C					(used when reading 1st record of file)
C				0 : adds one to record counter
C 	Output values of iProcessGen8_n:
C	iRec	integer		record # on file
C	cSrc	character	Source name
C	nYr	integer		year and ..
C	Doy	real		.. day of year of observation
C	XCsc	real		modified Carington variable for sub-SC point
C				  (NCoff already subtracted)
C	XCobs	real		modified Carrington variable of point P after
C				  trace-back to distance Radius (NCoff already subtracted)
C	xLat	real		heliographic latitude point P (deg)
C	xVV	real		IPS velocity (km/s)
C	xDVV	real		IPS velocity error (km/s)
C	xGG	real		IPS disturbance factor
C	xDGG	real		IPS disturbance factor error
C	xFRE	real		Source observation radio frequency
C	xSSS	real		IPS source size
C	dRA	real		RA(point P)-RA(Sun) (deg)
C	rLngSun real		geocentric ecliptic longitude Sun (deg)
C	rLngP	real		geocentric ecliptic longitude of point P (deg)
C					(relative to Sun)
C	rLatP	real		geocentric ecliptic latitude  of point P (deg)
C	rEloP	real		geocentric source elongation (deg)
C				NOTE: the sign is used to indicate East/West of Sun
C					rEloP>0: East of Sun; rEloP<0: West of Sun
C MODIFICATION HISTORY:
C	JUN-1994, Paul Hick (UCSD)
C	OCT-1996, Paul Hick (UCSD); added Ooty routines
C-
	subroutine ReadVIPSn8_n(iSys,iReadGen8_n,iProcessGen8_n,cWildVIPS,
     &		nCar,JDCar,bAuto,XCtest1,XCtest2,NCoff,MJDref,
     &		Radius,XElow,XEhigh,iEorW,XRlim,
     &		NLmax,NL,cvgfiles,SRCVG,SRC,SRCVGsav,SRCsav,
     &		IYRF,IREC,IYRS,DOYS8,XDS,XLS,XLL,XDL,XCE,XE,XC,YL,VV,DVV,GG,DGG,FRE,SSS,
     &		NSmax,IYRFsav,IRECsav,IYRSsav,DOYSsav8,XDSsav,XLSsav,XLLsav,
     &		XDLsav,XCEsav,XEsav,XCsav,YLsav,VVsav,DVVsav,GGsav,DGGsav,FRsav,SSSsav)

	real*8		JDCar(nCar),
     &			MJDref,
     &			xTT,JEpoch,
     &			Doy8,
     &			DOYS8(NLmax),DOYSsav8(NSmax)		! Time: doy of year

	integer		IYRF(NLmax), IYRFsav(NSmax),
     &			IREC(NLmax), IRECsav(NSmax),
     &			IYRS(NLmax), IYRSsav(NSmax)

	real		XDS (NLmax), XDSsav (NSmax),		! Sun-Earth distance
     &			XLS (NLmax), XLSsav (NSmax),		! Ecliptic longitude Sun
     &			XLL (NLmax), XLLsav (NSmax),
     &			XDL (NLmax), XDLsav (NSmax),
     &			XCE (NLmax), XCEsav (NSmax),
     &			XE  (NLmax), XEsav  (NSmax),		! Elongation point P
     &			XC  (NLmax), XCsav  (NSmax),		! Carrington variable of Point P
     &			YL  (NLmax), YLsav  (NSmax),		! Heliographic latitude of Point P
     &			VV  (NLmax), VVsav  (NSmax),		! Solar wind velocity (km/s)
     &			DVV (NLmax), DVVsav (NSmax),		! Solar wind velocity error (km/s)
     &			GG  (NLmax), GGsav  (NSmax),		! G-factor
     &			DGG (NLmax), DGGsav (NSmax),		! G-factor error
     &			FRE (NLmax), FRsav  (NSmax),		! Source observation frequency
     &			SSS (NLmax), SSSsav (NSmax)		! Source size

	character	cWildVIPS*100, cWildVIPSo*100,
     &			cvgfiles(*)*11,
     &			SRC(NLmax)*9,SRCsav(NSmax)*9,cSrc*9,SRCVG(NLmax)*115,SRCVGsav(NSmax)*115,cSRCVG*115,
     &			cSay*10			/'ReadVIPSn8'/,
     &			cFile*100,cFiles*11

	integer		iFrstYr			/0/,
     &			NS			/0/

	real		XCmargin		/1./

	logical		bAuto,
     &			bOpen			/.FALSE./,
     &			bGetRecord,
     &			bReading,
     &			bOpenFile,
     &			bMore,
     &			bReadErr,
     &			bSaveOverrun		/.FALSE./

	integer		nYrLen			/4/		! # digits used to identify year in file name
	integer		Str2Str, Flt2Str

	include		'dirspec.h'
	include		'filparts.h'
	include		'openfile.h'

	save		bOpen, bGetRecord, bSaveOverrun,
     &			iU, iYr, NS,
     &			iFrstYr, iFile

	logical		bLOSCheckLocn8_n
	integer		NLOS			/0/
	save		XCbeg,XCend,NLOS,NINSIDE,DLOS,cWildVIPSo,iRec1
C-------
C On first call check whether data files are present
C
	
	iRec1=0							! initialize the record length to 0 and save it.
	if (iSys.gt.0) then
	   print *,' '
	   if(abs(iSys).eq.3) print *, ' In readvipsn8_n. This is an Ooty IPS velocity read'
	   if(abs(iSys).eq.4) print *, ' In readvipsn8_n. This is a Nagoya IPS velocity read' 
	   if(abs(iSys).eq.5) print *, ' In readvipsn8_n. This is a KSWC IPS velocity read' 
	   if(abs(iSys).eq.6) print *, ' In readvipsn8_n. This is a MEXART IPS velocity read'
	   if(abs(iSys).eq.7) print *, ' In readvipsn8_n. This is a Pushchino IPS velocity read'
	   if(abs(iSys).eq.8) print *, ' In readvipsn8_n. This is a LOFAR IPS velocity read'
        else
	   print *,' '
	   if(abs(iSys).eq.3) print *, ' In readvipsn8_n. This is an Ooty IPS g-level read' 
	   if(abs(iSys).eq.4) print *, ' In readvipsn8_n. This is a Nagoya IPS g-level read' 
	   if(abs(iSys).eq.5) print *, ' In readvipsn8_n. This is a KSWC IPS g-level read' 
	   if(abs(iSys).eq.6) print *, ' In readvipsn8_n. This is a MEXART IPS g-level read'
	   if(abs(iSys).eq.7) print *, ' In readvipsn8_n. This is a Pushchino IPS g-level read'
	   if(abs(iSys).eq.8) print *, ' In readvipsn8_n. This is a LOFAR IPS g-level read'
	end if

	if(cWildVIPS.ne.cWildVIPSo) then
	  iFrstYr = 0	  
	  NS      = 0
	  cWildVIPSo = cWildVIPS
	end if
	if (iFrstYr .eq. 0) then
C	    print *, 'in readvipsn8_n', cWildVIPS

	    I = iSetFileSpec(cWildVIPS)
	    I = iGetFileSpec(FIL__NAME,FIL__TYPE,cFile)
	    iFile = index(cFile,cWildChar(:nYrLen))
	    if (iFile .eq. 0) call Say(cSay,'E','INVTMPLTE',
     &			'wildcard template must contain '//cWildChar(:nYrLen))

	    if (iSearch(1,cWildVIPS,cFile) .ne. 1)
     &		call Say(cSay,'E','FILNOTFOUND',
C     &			'no files matching wildcard found :'//cWildVIPS)
     &			'no files matching wildcard found :')

	    I = iSetFileSpec(cFile)
	    I = iGetFileSpec(FIL__NAME,FIL__TYPE,cFile)

	    I = 1

	      call Str2Flt_Int(.TRUE.)			! Save year of first file
	      call Str2Flt(cFile(iFile:iFile+nYrLen-1),I,iFrstYr)
	      iFile = index(cWildVIPS,cWildChar(:nYrLen)) ! Save position of wildcard %%%% in cWildVIPS

	end if
C------
C Check the save scratch arrays for valid data points.
C NOTE: NS can only be unequal to zero in AUTO mode. In AUTO mode only XCtest1
C	and XCtest2 change (both increase). All other data selection
C	criteria remain the same.

C	NL = 0				! Initialize data point counter Not needed to be set if set 0 outside (B. Jackson 10/04/2014)
	LOS = 0				! Counter for rejected lines of sight by bLOSCheckLocn8_n

	if (NS .ne. 0) then		! Only in AUTO mode
C-------
C Scan the `save' arrays. Only points with XC >= XCtest1 are needed for
C subsequent maps. Discard points with XC < XCtest1.

	    N = 0
	    do I=1,NS
		if (XCsav(I) .ge. XCtest1) then
		    N = N+1
		    SRCVGsav(N)= SRCVGsav(I)
		    SRCsav(N)  = SRCsav(I)
		    IYRFsav(N) = IYRFsav(I)
		    IRECsav(N) = IRECsav(I)

		    IYRSsav(N) = IYRSsav(I)
		    DOYSsav8(N)= DOYSsav8(I)
 		    XDSsav (N) = XDSsav (I)
		    XLSsav (N) = XLSsav (I)
		    XLLsav (N) = XLLsav (I)
		    XDLsav (N) = XDLsav (I)
		    XCEsav (N) = XCEsav (I)
		    XEsav  (N) = XEsav  (I)
		    XCsav  (N) = XCsav  (I)
		    YLsav  (N) = YLsav  (I)
		    VVsav  (N) = VVsav  (I)
		    DVVsav (N) = DVVsav (I)
		    GGsav  (N) = GGsav  (I)
		    DGGsav (N) = DGGSav (I)
		    FRsav  (N) = FRsav  (I)
		    SSSsav (N) = SSSSav (I)
		end if
	    end do
	    NS = N			! Valid # point in save arrays

C--------
C Pick up all points in the save arrays which are needed for the current
C map. Since all points with XC < XCtest1 have just been discarded, the
C points are only tested for XC <= XCtest2.

	    do I=1,NS
		if (XCsav(I) .le. XCtest2) then
		    NL = NL+1
	            IL1 = iSetFileSpec(cFile)
	            IL1 = iGetFileSpec(FIL__NAME,FIL__TYPE,cFiles)
	    	    cvgfiles(NL) = cFiles
		    SRCVG(NL)= SRCVGsav(I)
		    SRC(NL)  = SRCsav(I)
		    IYRF(NL) = IYRFsav(I)
		    IREC(NL) = IRECsav(I)

		    IYRS(NL) = IYRSsav(I)
		    DOYS8(NL)= DOYSsav8(I)
 		    XDS (NL) = XDSsav (I)
 		    XLS (NL) = XLSsav (I)
 		    XLL (NL) = XLLsav (I)
 		    XDL (NL) = XDLsav (I)
 		    XCE (NL) = XCEsav (I)
 		    XE  (NL) = XEsav  (I)

		    XC  (NL) = XCsav  (I)
		    YL  (NL) = YLsav  (I)
		    VV  (NL) = VVsav  (I)
		    DVV (NL) = DVVsav (I)
		    GG  (NL) = GGsav  (I)
		    DGG (NL) = DGGsav (I)
		    FRE	(NL) = FRsav  (I)
		    SSS	(NL) = SSSsav (I)
		endif
	    end do
	end if

C-------
C If no file open then check which file is to be read.
C (bOpen=.TRUE. only happens in AUTO mode. Then also bGetRecord=.FALSE.,
C indicating an unprocessed record is still in memory).
C The data files are assumed to be organized in yearly files: one file
C for each year of data. XCtest1-XCmargin is used to determine which yearly
C file must be opened first.

	if (.not. bOpen) then

	    XCt = XCtest1-XCmargin
	    I = XCt
	    if (I .lt. 1 .or. I .ge. nCar) call Say(cSay,'E','InvXCtest1',
     &		    'start value outside range covered by JDCar')

	    xTT = JDCar(I)+(XCt-I)*(JDCar(I+1)-JDCar(I))
	    call Julian8(1,iYr,Doy8,xTT,JEpoch)
	    iYr = max(iYr,iFrstYr)	! Year of file to be opened

	    iYr = iYr-1			! Will be incremented again in a second
	    bGetRecord = .TRUE.		! Make sure to read record

	end if
C-------
C  Read data until
C  -  end point XCtest2 is reached, or
C  -  EOF is reached, or
C  -  arrays are full
C
C bOpen and bGetRecord are used in subsequent calls.
C If a record has been read succesfully, but is rejected because it's out
C of range (iProc=2) or because the arrays are full (NL>=NLmax), then
C the file is kept open (bOpen=.TRUE.) and bGetRecord is set .FALSE.
C On the next call to ReadVIPSn8_n, the procedure will skip straight to the
C iProcessGen8_n call to try processing the last record read in the previous
C call before reading the next record from the file that is still open.

	bMore	 = .TRUE.
	bReading = .FALSE.		! Only used to control messages to screen

	do while (bMore)

	    if (.not. bOpen) then	! Open next year's file
		iYr = iYr+1
		iFirst = 1		! Force ReadVIPSn8_n to initialize record counter
C-------
C Complete the data file name by changing the file extension to the
C proper year. Open data file. If bOpenFile returns 0 the file probably
C didn't exist (exit do while loop). bOpenFile assigns iU if successfully
C opening the file.

		cFile = cWildVIPS	! Replace wildcard by year
		I = Int2Str(iYr,cFile(iFile:iFile+nYrLen-1))
		
		bOpen = bOpenFile(OPN__REOPEN+OPN__TEXT+OPN__READONLY,iU,cFile,iRecl)
		bMore = bOpen		! Open error will exit while loop
	        IL1 = iSetFileSpec(cFile)
	        IL1 = iGetFileSpec(FIL__NAME,FIL__TYPE,cFiles)
	    end if

	    bReadErr = .FALSE.
	    if (bOpen .and. bGetRecord) then
		if (.not. bReading) then
		    call Say(cSay,'I','Reading',cFile)
		    bReading = .TRUE.
		end if
		bMore = iReadGen8_n(iSys,iU,iFirst) .eq. 0

		if (bMore) then
		    iFirst = 0		! Record counter has been initialized
		else			! Read error: EOF assumed
		    iU = iFreeLun(iU)
		    bOpen = .FALSE.
		    bReadErr = .TRUE.
		end if

	    end if
	    if (bReadErr) then		! After read error, open new file
		bMore = .TRUE.
	    else if (bMore) then

		iProc = iProcessGen8_n(iR,nCar,JDCar,NCoff,XCtest2,Radius,
     &			cSRCVG,cSrc,nYr,Doy8,XCsc,XCobs,xLat,xVV,xDVV,xGG,xDGG,xFRE,xSSS,dRA,
     &			rLngSun,rDisSun,rLngP,rLatP,rEloP)

		bMore = iProc .ne. 2		! Check for XCsc past XCtest2
		if (iProc .eq. 1) then		! Record contains valid data with XCsc <= XCtest2

		    call Julian8(10,nYr,Doy8,xTT,JEpoch)
						! All XE (East+West), only positive XE (East) or only negative XE (West)
		    Rdum = (1-abs(iEorW))*abs(rEloP)+iEorW*rEloP

     		    if ((XElow .le. Rdum .and. Rdum .le. XEhigh) .and.
     &			( abs(dRA) .le. XRlim ) .and.
     &			(MJDref .eq. 0 .or. xTT .le. MJDref)
     &		       ) then

			if (XCtest1 .le. XCobs) then

			    if (bAuto) then
			        bSaveOverrun = NS .ge. NSmax
				if (.not. bSaveOverrun) then
				    NS = NS+1
				    SRCVGsav(NS)= cSRCVG
				    SRCsav(NS)	= cSrc
				    IYRFsav(NS) = iYr
				    IRECsav(NS) = iR

				    IYRSsav(NS) = nYr
				    DOYSsav8(NS)= Doy8
				    XDSsav (NS) = rDisSun
				    XLSsav (NS) = rLngSun
				    XLLsav (NS) = rLngP
				    XDLsav (NS) = rLatP
				    XCEsav (NS) = XCsc
				    XEsav  (NS) = rEloP

				    XCsav  (NS) = XCobs
				    YLsav  (NS) = xLat
				    VVsav  (NS) = xVV
				    DVVsav (NS) = xDVV
				    GGsav  (NS) = xGG
				    DGGsav (NS) = xDGG
				    FRsav  (NS) = xFRE
				    SSSsav (NS) = xSSS
				end if

			    end if
C--------
C Pick up all points needed for the current map.Test only for XC <= XCtest2
			    if (XCobs .le. XCtest2) then
C--------
C Check what part of line of sight is inside the actual boundaries [XCbeg,XCend]

				if (bLOSCheckLocn8_n(XCbeg,XCend,NLOS,NINSIDE,DLOS,
     &						nYr,Doy8,rLngP,rLatP,XCsc,rLngSun,rDisSun)) then

				    if (NL .lt. NLmax) then

					NL = NL+1
	    	    			cvgfiles(NL) = cFiles
				        SRCVG(NL)= cSRCVG
					SRC(NL)	 = cSrc
					IYRF(NL) = iYr
					IREC(NL) = iR

					IYRS(NL) = nYr		! Time of observation
					DOYS8(NL)= Doy8
					XDS (NL) = rDisSun	! Sun-Earth distance (AU)
					XLS (NL) = rLngSun	! Geocentric ecliptic longitude of Sun (deg)
					XLL (NL) = rLngP	! Geocentric ecliptic lng(P)-lng(Sun)
					XDL (NL) = rLatP	! Geocentric ecliptic lat(P)
					XCE (NL) = XCsc		! Carrington variable sub-Earth point
					XE  (NL) = rEloP
					XC  (NL) = XCobs

					YL  (NL) = xLat
					VV  (NL) = xVV
					DVV (NL) = xDVV
					GG  (NL) = xGG
					DGG (NL) = xDGG
					FRE (NL) = xFRE
					SSS (NL) = xSSS

C	write (*,'(I4,I5,9F10.3)') NL,IYRS(NL),Doy8,XC(NL),YL(NL),VV(NL),DVV(NL),GG(NL),DGG(NL),FRE(NL),SSS(NL)
C	write (*,'(I4,F6.2,4F10.2)') ID,UT,XCtest1,XCobs,XCtest2

				    else
					call Say(cSay,'I','NLmax','Data arrays full: not all data processed')
					bMore = .FALSE.
				    end if
				else
				    LOS = LOS+1

				end if
			    end if
			end if
		    end if
		end if

C-------
C At this point bMore = .FALSE. if iProc=2 (XCsc past XCtest2) or the data arrays are full

		bGetRecord = bMore .or. .not. bGetRecord

	    end if

	end do

	if (.not. bAuto .and. bOpen) then	! If not in AUTO mode, make sure to close open file
	    iU = iFreeLun(iU)
	    bOpen = .FALSE.
	end if

	if (bSaveOverrun) then			! Can only happen in AUTO mode
	    call Say(cSay,'W','NSmax','Save arrays overrun: contents discarded')
	    NS = 0
	    bSaveOverrun = .FALSE.
	end if

	I = 0
	I = I+Int2Str(NL,cFile(I+1:))+1
	I = I+Str2Str('sources in',cFile(I+1:))+1
	I = I+Flt2Str(NCoff+XCtest1,2,cFile(I+1:))
	I = I+Str2Str('-',cFile(I+1:))
	I = I+Flt2Str(NCoff+XCtest2,2,cFile(I+1:))

	if (NLOS .ne. 0) then
	    I = I+2
	    I = I+Int2Str(LOS,cFile(I+1:))+1
	    I = I+Str2Str('; sources rejected',cFile(I+1:))
	end if
	call Say(cSay,'I','# obs',cFile)

	return

	entry ReadVipsLOSCheckn_n(XCbegIn,XCendIn,NLOSIn,NINSIDEIn,DLOSIn)

	XCbeg 	= XCbegIn
	XCend 	= XCendIn
	NLOS  	= NLOSIn
	NINSIDE = NINSIDEIn
	DLOS  	= DLOSIn

	return
	end

	logical function bLOSCheckLocn8_n(XCbeg,XCend,NLOS,NINSIDE,DLOS,iYr,Doy8,rLngLOS,rLatLOS,ObsXC,rLngSun,rDisSun)

	real*8		Doy8
	
	bLOSCheckLocn8_n = NLOS .eq. 0
	if (bLOSCheckLocn8_n) return

	ObsLng = rLngSun-180.

	Inside = 0
	do J=1,NLOS			! Loop over all segments along LOS
	    P1 = rLngLOS		! Topocentric ecliptic longitude diff. with Earth-Sun line
	    P2 = rLatLOS		! Topocentric ecliptic latitude
	    P3 = (J-.5)*DLOS		! Distance along LOS (units of ObsR)
					! Convert to heliocentric coordinates
	    call POINT_ON_LOS(P1,P2,P3,ELO,ELOS,iEorW)
	    P3 = P3*rDisSun		! Convert to AU

	    P1 = ObsLng+P1		! Heliocentric longitude of P
	    call ECLIPTIC_HELIOGRAPHIC8(0,iYr,Doy8,P1,P2)
					! Convert to heliographic coordinates
	    P1 = XMAP_OBS_POS(ObsXC,P1)	! Modified Carrington variable of P

	    if (XCbeg .le. P1 .and. P1 .le. XCend) Inside = Inside+1
	end do

	bLOSCheckLocn8_n = Inside .ge. NINSIDE

	return
	end
C+
C NAME:
C	iReadGen8_n
C PROCEDURE:
C >	Structure of Gen files from Nagoya:
C	iY,iM,iD	integer*8 	year, month and day of observation in UT (yyyymmdd)
C	ih,im,is	integers	hour, minute, sec of observation start in UT (hh:mm:ss)
C	Site		string		STEL327
C	cSrce		string		descriptor for radio source
C	iRA(3)		integers	hour, minute, sec in J2000 coordinate
C	  save		iY,iMO,iD,iH,iMI,UT,cSrce,iRA50,iDEC50,iRA,
C     &			iDEC,ELO,FREQ,SSIZE,iVEL,SN,SCINDX,GVAL,iRec1,cSRCVGe
C	  save		DP,iHLA,iHLO,iGLA,iGLO,iROT,iDVEL

C              cFmt*40 '(A9,3I2,F6.2,F5.2,I4,I5,I4,3I5,I4,F8.5)'/s (hh mm ss)
C	iDEC(3)		integers	deg, minute, sec in J2000 coordinates (deg mm ss)
C	Limb		string		solar limb of source (E - east, W - west)
C	Dist		real		P-point distance (in solar radii)
C	Lat.		real		Latitude of radio source p-point (in degrees)
C	PA		real		Ecliptic Position angle  (in degrees)
C	Elong		real		elongation of the radio source (in degrees)
C	iVel1		integer		IPS single-site velocity (km/s)
C	idVel1		integer		IPS single-site velocity error (km/s)
C	gval1		real		IPS single-site g-value (at Toyokawa)
C	gdval1		real		IPS single-site g-value error (at Toyokawa)
C	Method1		string		descriptor indicating single or multi site
C	iVel3		integer		IPS three-site velo	read (iU,'(4I,2I2,1X,3(I2,1X)city (km/s)
C	idVel3		integer		IPS single-site velocity error (km/s)
C	gval3		real		IPS single-site g-value (at Toyokawa)
C	gdval3		real		IPS single-site g-value error (at Toyokawa)
C	Method3		string		descriptor indicating single or multi site
C
C These are set up currently for the new Nagoya files. When we have a standard format, hopefully Ooty, MEXART, and the KSWC will follow this. 
C-
	function iReadGen8_n(iSys,iU,iFirst)
	character	cMon*3,
     &			cSrce*8,cSrc*9,cSRCVGe*115,cSRCVG*115
	character	cInput1*166
	character	onespace*1/' '/
	character	name1*4
	character*2	cshortcite
	character*3	cRA(3),cDEC(3)
	character*4	ciVel1						! This may be changed to an integer some day
	character*5	cSite
	character*7	cdiVel1,cgdval1,cdgval3				! These may all be changed to integers some day
	real*8		JDCar(nCar),
     &			dLngSun, dLatSun, dDisSun
     	real*8		Doy8
	real*8		JD1,JEpoch1,JD2,JEpoch2,Doy18,Doy28
	integer		iRA (3), iRA2(3),
     &			iDEC(3),iDEC2(3)
	external	EARTH8
	save		iY,iMO,iD,iH,iMI,UT,cSrce,iRA50,iDEC50,iRA,
     &			iDEC,ELO,FREQ,SSIZE,iVEL,iDVEL,SN,SCINDX,GVAL,dGVAL,iRec1,iSYS1,cSRCVGe
	save		DP,iHLA,iHLO,iGLA,iGLO,iROT
	include		'sun.h'

C 	New general format version

	read (iU,'(A166)',iostat=iReadGen8_n) cInput1					! read CELIAS file headers	
	read (cInput1,'(A4)') name1

	if(name1.ne.'Date') then
C	  print *, 'beginng read of g-level data'
	  read (cInput1(1:166),'(I4,2I2,1X,3(I2,1X),5X,A5,F4.0,4X,A8,I5,1X,3I3,1X,3I3,32X,2F7.0,2F7.3,9X,2F7.0,2F7.3)'),
     &		iY,iMO,iD,iH,iMI,iSI,cSite,Freq,cSrce,iSize,iRA2(1),iRA2(2),iRA2(3),iDEC2(1),iDEC2(2),IDEC2(3),
     &		vel3,dvel3,gval3,gdval3,vel1,dvel1,gval1,dgval1

C	  print *, 'Now writing the general formatted data'
C	  write(*,'(I4,2I2,1X,3(I2,1X),5X,A5,F4.0,4X,A8,I5,1X,3I3,1X,3I3,32X,2F7.0,2F7.3,9X,2F7.0,2F7.3)'),
C     &		iY,iMO,iD,iH,iMI,iSI,csite,Freq,cSrce,iSize,iRA2(1),iRA2(2),iRA2(3),iDEC2(1),iDEC2(2),IDEC2(3),
C     &		vel3,dvel3,gval3,gdval3,vel1,dvel1,gval1,dgval1
C
C	  Precess from IRA2000, DEC2000 to present, end up with the present-day IRA, IDEC
C
	  UT  = float(iH) + float(iMI)/60.0 + float(iSI)/3600.0				! This is the midpoint of the observation, UT in Hr
 	  iY2 = iY		
	  call Precession_Drive(3,2000,0.0D0,iY2,0.0D0,iMO,iD,UT,iRA2,iDEC2,iRA,iDEC,cRA,cDEC)
C
C         if(iSys.eq.3.or.iSys.eq.4.or.iSys.eq.5.or.iSys.eq.6.or.iSys.eq.7.or.iSys.eq.8) then  
C             ! Ooty, Nagoya, KSWC, MEXART, Pushchino or LOFAR Data - standard format

	  if(iSys .eq. 3 .or. iSys .eq. -3) then 								! Ooty standard format 
C	    print *, 'This is an Ooty file in the new format' 
	    iVEL  = nint(vel3)
	    iDVEL = nint(dvel3)
	    GVAL  = gval1
	    dGVAL = dgval1
	    if (nint(Freq).ne.327) then
	      write(*,'(A,I6,I5,2I3,F7.3,E9.3,I6)') 'Radio frequency for Ooty source bad =', iRec1,iY2,iMO,iD,UT,Freq,nint(Freq)
	      Freq = 327.0
	    end if
	    cshortsite = 'O '
	    if(iSize.eq.-999) SSIZE = 0.10
	  end if

	  if(iSys .eq. 4 .or. iSys .eq. -4) then 								! Nagoya standard format 
C	    print *, 'This is a Nagoya file in the new format' 
	    iVEL  = nint(vel3)
	    iDVEL = nint(dvel3)
	    GVAL  = gval1
	    dGVAL = dgval1
	    if (nint(Freq).ne.327) then
	      write(*,'(A,I6,I5,2I3,F7.3,E9.3,I6)') 'Radio frequency for STELab source bad =', iRec1,iY2,iMO,iD,UT,Freq,nint(Freq)
	      Freq = 327.0
	      stop
	    end if
	    cshortsite = 'S '
	    if(iSize.eq.-999) SSIZE = 0.10
	  end if

	  if(iSys .eq. 5 .or. iSys .eq. -5) then 								! Nagoya standard format 
C	    print *, 'This is a KSWC file in the new format' 
	    iVEL  = nint(vel3)
	    iDVEL = nint(dvel3)
	    GVAL  = gval1
	    dGVAL = dgval1
	    if (nint(Freq).ne.327) then
	      write(*,'(A,I6,I5,2I3,F7.3,E9.3,I6)') 'Radio frequency for KSWC source bad =', iRec1,iY2,iMO,iD,UT,Freq,nint(Freq)
	      Freq = 327.0
	      stop
	    end if
	    cshortsite = 'K '
	    if(iSize.eq.-999) SSIZE = 0.10
	  end if

	  if(iSys .eq. 6 .or. iSys .eq. -6) then 								! MEXART standard format 
C	    print *, 'This is a MEXART file in the new format' 
	    iVEL  = nint(vel1)
	    iDVEL = nint(dvel1)
	    GVAL  = gval1
	    dGVAL = dgval1
	    if (nint(Freq).ne.140) then
	      write(*,'(A,I6,I5,2I3,F7.3,E9.3,I6)') 'Radio frequency for MEXART source bad =', iRec1,iY2,iMO,iD,UT,Freq,nint(Freq)
	      stop
	    end if
	    cshortsite = 'M '
	    if(iSize.eq.-999) SSIZE = 0.20
	  end if

	  if(iSys .eq. 7 .or. iSys .eq. -7) then 								! Pushchino standard format 
C	    print *, 'This is a Pushchino file in the new format' 
	    iVEL  = nint(vel1)
	    iDVEL = nint(dvel1)
	    GVAL  = gval1
	    dGVAL = dgval1
	    if (nint(Freq).ne.110) then
	      write(*,'(A,I5,I5,2I3,F7.3,E9.3,I6)') 'Radio frequency for Pushchino source bad =', iRec1,iY2,iMO,iD,UT,Freq,nint(Freq)
	      Freq = 110.0
C	      stop
	    end if
	    cshortsite = 'P '
	    if(iSize.eq.-999) SSIZE = 0.25
	  end if

	  if(iSys .eq. 8 .or. iSys .eq. -8) then 								! LOFAR standard format 
C	    print *, 'This is a LOFAR file in the new format' 
	    iVEL  = nint(vel3)
	    iDVEL = nint(dvel3)
	    GVAL  = gval1
	    dGVAL = dgval1
	    cshortsite = 'L '
	    if(iSize.eq.-999) SSIZE = 0.20
	  end if

	  if(iy.lt.2000) then
	     iYRM = iY - 1900
	  else
	     iYRM = iY - 2000
	  end if

C	  write(*,'(A8,1X,3I2,F6.2,1X,3A3,1X,3A3,1X,A2,F5.3,I5,I5,2F7.3)')
C     &		cSrce,iYRM,iMO,iD,UT,cRA(1),cRA(2),cRA(3),cDEC(1),cDEC(2),cDEC(3),cshortsite,SSIZE,iVEL,iDVEL,GVAL,dGVAL

	  write(cSRCVGe(1:73),'(A8,1X,3I2,F6.2,1X,3A3,1X,3A3,1X,A2,F5.3,I5,I4,F8.5,F7.3)')
     &		cSrce,iYRM,iMO,iD,UT,cRA(1),cRA(2),cRA(3),cDEC(1),cDEC(2),cDEC(3),cshortsite,SSIZE,iVEL,iDVEL,GVAL,dGVAL

	  iRec1 = (1-iFirst)*iRec1+1		! iFirst=1 initializes record counter
	end if
	if(iSys1.ne.iSys) iRec1 = 0
	iSys1 = iSys
	return
C+
C NAME:
C	iProcessGen8_n
C CALLING SEQUENCE:
	entry iProcessGen8_n(iR,nCar,JDCar,NCoff,XCend,Radius,
     &		cSRCVG,cSrc,nYr,Doy8,XCsc,XCobs,xLat,xVV,xDVV,xGG,xDGG,xFRE,xSSS,dRA,
     &		rLngSun,rDisSun,rLngP,rLatP,rEloP)
C-
C	print *, 'iVEL =',iVEL
C	if (iVEL .eq. 0) then			! No velocity available: reject
C	    iProcessGen8_n = 0		
C	    return				! We can't do this check if the same entry is used for density - BVJ 9/16/2014
C	end if

	if (irec1 .eq. 0) then			! No record length set: reject
	    iProcessGen8_n = 0		
	    return
	end if

	cSRCVG = cSRCVGe
	cSrc	= '         '
	cSrc = cSrce
	nYr  = iY

	cMon = '   '				! Make sure iMO is used
	call DATE_DOY(0,nYr,cMon,iMO,iD,iDoy)
	Doy8 = iDoy+dble(UT/24.)

	XCsc = XMAP_SC_POS8(EARTH8,nYr,Doy8,-nCar,JDCar)

	if (XCsc .eq. BadR4()) then		! Time outside JDCar range
		iProcessGen8_n = 0
		return
	end if
						! Sub-Earth position
	if (XCsc .gt. XCend) then		! Sub-Earth past XCend
	    iProcessGen8_n = 2
	    irec1 = 0				! Needed in case the first line of the file is a header
	    return
	end if

	iR   = iRec1

	call SunNewcomb8(0,nYr,Doy8,dLngSun,dLatSun,dDisSun)
						! Geocentric eclip. coord(Sun)
	rDisSun = dDisSun
	rLngSun = dLngSun			! Geocentric eclip. lng(Sun)

	rLng = dLngSun
	rLat = dLatSun				! Geocentric ecliptic lat(Sun)
	call ECLIPTIC_EQUATOR8(0,nYr,Doy8,rLng,rLat)
	dRA = rLng				! RA (Sun)

	rLng = (float(iRA(1))+(float(iRA(2))+float(iRA(3))/60.)/60.)*15.	! RA (P) (deg)
	if(iDEC(1).lt.0.or.iDEC(2).lt.0.or.iDEC(3).lt.0) then
	  if(iDEC(3).lt.0) then
	    rLat = float(iDEC(3))/3600.0
	  else
	    if(iDEC(2).lt.0) then
	      rLat = float(iDEC(2))/60.0 - float(iDEC(3))/3600.0
	    else
	      rLat = float(iDEC(1)) - float(iDEC(2))/60.0 - float(iDEC(3))/3600.0
	    end if
	  end if
	else
	  rLat = float(iDEC(1)) + (float(iDEC(2)) + float(iDEC(3))/60.)/60.	! Decl (P) (deg)
	end if

C	rLng  = (iRA(1)+(iRA(2)+iRA(3)/60.)/60.)*15.	! RA (P) (deg)
C	rLat  = iDEC(1)+(iDEC(2)+iDEC(3)/60.)/60.	! Decl (P)

	dRA = rLng-dRA				! RA(P)-RA(Sun)
	if (abs(dRA) .gt. 180.) dRA = dRA-sign(1.,dRA)*360.

	call ECLIPTIC_EQUATOR8(1,nYr,Doy8,rLng,rLat)
						! Geocentric eclip. coord(P)
	rLng  = rLng-rLngSun
	rLngP = rLng				! Geocentric ecliptic lng(P)-lng(Sun)
	rLatP = rLat				! Geocentric ecliptic lat(P)

	rP = max(0.25, cosd(rLng)*cosd(rLat) )	! Geocentric distance P (at least 0.25 AU)
	call POINT_ON_LOS(rLng,rLat,rP,rEloP,EloSun,i)
						! Heliocentric coord(P) rel. Sun
C	print *, ELO,rEloP, ELO-rEloP

	rP    = rP*dDisSun			! Dist. Sun-P (AU)
	rEloP = i*rEloP				! rEloP>0 : East of Sun
						!      <0 : West of Sun
	rLng = rLngSun+180.+rLng		! Heliocentric eclip. lng(P)
	call ECLIPTIC_HELIOGRAPHIC8(0,nYr,Doy8,rLng,rLat)
						! Heliographic coord(P)
	xLat = rLat				! Heliographic lat(P)

	XCobs = XMAP_OBS_POS(XCsc,rLng)		! Point-P at rP
	XCobs = XCobs+(Radius-rP)/iVEL*SUN__SPIRAL ! Point-P traced to Radius

	xVV  = iVEL				! IPS velocity
	xDVV = iDVEL				! IPS velocity error
	xGG  = GVAL				! IPS g-value
	xDGG = dGVAL				! IPS g-value error
	xFRE = FREQ				! IPS frequency for this source
	xSSS = SSIZE				! IPS source size

	iProcessGen8_n = 1

	return
	end
