C+
C NAME:
C	Average
C PURPOSE:
C	Replaces groups of points with their averages if those points are 
C	closer together in time than a specified value.
C CALLING SEQUENCE:
	program Average
C INPUTS:
C	Helios data file
C	The data for each photometer/color/filter combinations must be
C	chronologically ordered
C OUTPUTS:
C	Averaged Helios data file
C CALLS:
C	BadR4, iHOSInfo, Say, AskR4, AskChar, AskYN, bOpenFile, iFreeLun, iHOSRead, iHOSWrite
C	iPutFileSpec, iGetFileSpec, ArrR4Zero, ArrI4Zero, ArrR4DivideByArrI4
C INCLUDE:
	include		'hos_e9.h'
	include		'filparts.h'
	include		'openfile.h'
C PROCEDURE:
C	The input file does not need to be rigourously chronological, but should be
C	chronologically ordered for each individual photometer/color/filter combination.
C
C	The output file also will not be rigourously chronological. Individual 
C	photometer/color/filter combination are stored in chronological order, and
C	are written to file sequentially. Combinations are processed in the triple
C	loop with the filter in the innermost loop, then the photometer, and the color
C	in the outermost loop.
C
C	If the 'merge colors' option is selected then the color values are ignored
C	(essentially the file is treated as if it contains only visual (color=3) data.
C	In the output file all color entries are set to 3.
C MODIFICATION HISTORY:
C	MAY-2000, Paul Hick (UCSD/CASS); pphick@ucsd.edu), rewritten and documented
C-
	parameter	(nS = 32)

	integer		Lst(3,11)

	integer		P
	integer		C
	integer		F
	integer		PP
	integer		CC
	integer		FF

	integer		CBeg
	integer		CEnd

	integer		kPP(3)		/HOS__P_1,HOS__P_2,HOS__P_3/
	integer		kCC(3)		/HOS__C_1,HOS__C_2,HOS__C_3/
	integer		kFF(5)		/HOS__F_1,HOS__F_2,HOS__F_3,HOS__F_4,HOS__F_5/

	integer		NZZ(nS)
	integer		FBegP(3)	/1,1,4/		! Start filters for all three photometers
	integer		FEndP(3)	/5,5,4/		! End filters

	real		SZZ(nS)
	real		ZZ (nS)
	real		dT		/0.0417/	! Default is one hour
	real		LL

	character	cSay*7		/'Average' /

	logical		bAllColors
	logical		bOpenFile

	character	cFileIn *80
	character	cFileOut*80

	Bad = BadR4()

	if (iHOSInfo(OPN__STOP,cFileIn,iRecl,iSc,iYr,TMin,TMax,kPCFAll,Lst) .eq. 0)
     &		call Say(cSay,'I','StopPlay','Exit')

	iS = 2					! Phot 3 'sectors': intensity and pB
	if (iRecl .gt. 6) iS = iRecl-4.5	! Sectors in file

	N = iPutFileSpec(FIL__TYPE,FIL__TYPE,'.ave')
	N = iGetFileSpec(0,FIL__TYPE,cFileOut)
	call AskChar('Output file$parse',cFileOut)

	call AskR4('Averaging interval (in days)',dT)
	call AskYN('Combine all colors$no',bAllColors)

	if (bOpenFile(OPN__HOS+OPN__READONLY+OPN__REOPEN+OPN__STOP,iUIn ,cFileIn ,iRecl)) continue
	if (bOpenFile(OPN__HOS+OPN__NEW     +OPN__REOPEN+OPN__STOP,iUOut,cFileOut,iRecl)) continue

	CBeg = 1
	if (bAllColors) CBeg = 3
	CEnd = 3

	TTOld = 0		! Better to initialize here. The value doesn't matter.

	kPCF = kPCFAll-iand(kPCFAll,HOS__PCF_ALL)

	iW = 0
	do C=CBeg,CEnd
	    do P=1,3
		do F=FBegP(P),FEndP(P)

		    kPCFIn = kPCF+kPP(P)+kFF(F)
		    if (bAllColors) then
			kPCFIn  = kPCFIn+HOS__C_ALL
		    else
			kPCFIn  = kPCFIn+kCC(C)
		    end if

		    N  = 0
		    iR = 0
		    do while (iHOSRead(kPCFIn,iUIn,iRecl,iR,TT,PP,CC,FF,RR,LL,iS,ZZ) .eq. 0)

			!-------
			! First check whether an averaged record needs to written to the output file.
			! This is the case if the current time TT is more than dT later than the previous time TTOld.
			! (N = 0 happens only for the 1st record, i.e. when entering the while loop).

			if (N .ne. 0 .and. TT-TTOld .gt. dT) then
			    STT = STT/N
			    SRR = SRR/N
			    SLL = SLL/N
			    call ArrR4DivideByArrI4(iS,SZZ,NZZ,SZZ)

			    if (bAllColors) CC = CBeg
			    if (iHOSWrite(kPCFAll,iUOut,iRecl,iW,STT,PP,CC,FF,SRR,SLL,iS,ZZ) .ne. 0)
     &				call Say(cSay,'E',cFileOut,'error writing to file')

			    N = 0				! Start collecting data for the next averaged record
			end if

			!-------
			! For N = 0 initialize variables for summing. N is zero when the very
			! 1st record is written to file, and each time a record is written to file.

			if (N .eq. 0) then			! Initialize
			    STT = 0
			    SRR = 0
			    SLL = 0
			    call ArrI4Zero(iS, NZZ)
			    call ArrR4Zero(iS, SZZ)
			end if

			!-------
			! Always sum the current record when N=0.
			! Otherwise sum only if the current time TT is closer than
			! dT to the previous time TTOld

			if (N .eq. 0 .or. TT-TTOld .le. dT) then
			    N   = N  +1
			    STT = STT+TT
			    SRR = SRR+RR
			    SLL = SLL+LL
			    do I=1,iS
				if (ZZ(I) .ne. Bad) then
				    NZZ(I) = NZZ(I)+1
				    SZZ(I) = SZZ(I)+ZZ(I)
				end if
			    end do
			end if

			TTOld = TT				! Current time becomes previous time

		    end do

		    if (N .ne. 0) then				! Write residual data
			STT = STT/N
			SRR = SRR/N
			SLL = SLL/N
			call ArrR4DivideByArrI4(iS,SZZ,NZZ,SZZ)

			if (bAllColors) CC = CBeg
			if (iHOSWrite(kPCFAll,iUOut,iRecl,iW,STT,PP,CC,FF,SRR,SLL,iS,ZZ) .ne. 0)
     &				call Say(cSay,'E',cFileOut,'error writing to file')
		    end if

		end do
	    end do
	end do

	iUIn  = iFreeLun(iUIn )
	iOOut = iFreeLun(iUOut)

	N = iHOSInfo(OPN__REOPEN,cFileOut,iRecl,iSc,iYr,TMin,TMax,kPCFAll,Lst)

	call Say(cSay,'I','StopPlay','Done')

	end