program dailym
c
c	a program calc daily means grib -> grib
c
c
c	file 50: seq grib file-x
c	file 51: index grib file-x
c	file 54: list of times [yymmdd00]
c	file 55: var <cr> level
c
c	file 1: output file seq grib
c
c	v 1.1 7/29/94 minor mod to grib desc .. added nsum to ave of anal/fcst
c		change to setbit(....) .. assumes mask in constant in time
c	v 1.2 3/95 assumes mask can change in time needed for cloud temp
c
        parameter (maxv=100000)

        logical lbms(maxv), hasbit
        real data(maxv), mean(maxv)
        integer icount, count(maxv)

        integer kpds(100),kgds(100),recl,offset
        integer kpds6, kpds7

        integer id, time, time2, year, mon, hour, day
        real level
        character name*8

        read(55,'(a8)') name
        read(55,*) level
        close(55)

c	find kpds(5..7) parameters

        call atok5(name,id)
        if (id.eq.-1) then
            write(*,*) 'illegal variable name ', name
            call exit(8)
        endif

        call setlev(kpds,kgds,level)
        kpds6=kpds(6)
        kpds7=kpds(7)

	write(*,*) 'dailym ',name,'->',id,' level=',level

        inum = 0
        icount = 0

250	read(54,'(i6)',end=500) time
            if (time.le.0) goto 500

            time2 = time
            day = mod(time,100)
            time = time/100
            mon = mod(time,100)
            time = time/100
            year = mod(time,100)

            do 255 i = 1, maxv
		count(i) = 0
                mean(i) = 0.0
255         continue
            nsum = 0

            do 480 hour = 0, 18, 6

                do 260 i = 1, 100
                    kpds(i) = -1
260             continue
                kpds(5) = id
                kpds(6) = kpds6
                kpds(7) = kpds7
                kpds(8) = year
                kpds(9) = mon
                kpds(10) = day
                kpds(11) = hour

                call pdscan(51,kpds,inum,offset,recl)
                if (recl.eq.0) then
                    write(*,*) 'missing time ', time2, hour
		    call exit(8)
                    goto 480
                endif

                call rdgb(50,recl,offset,kpds,kgds,ndata,lbms,data)

                call getbit(kpds, kgds, hasbit)

                nx = kgds(2)
                ny = kgds(3)
                nxny = nx * ny

                if (hasbit) then
                    do 270 i = 1, nxny
                        if (lbms(i)) then
			    mean(i) = mean(i) + data(i)
			    count(i) = count(i) + 1
			endif
270                 continue
                else
                    do 280 i = 1, nxny
                        mean(i) = mean(i) + data(i)
			count(i) = count(i) + 1
280                 continue
                endif
                nsum = nsum + 1
480         continue

            if (nsum.eq.4) then
c               *** find daily ave and write grib record ***
                do 490 i = 1, nxny
		    hasbit = .false.
		    if (count(i).gt.0) then
			mean(i) = mean(i) / float(count(i))
			lbms(i) = .true.
		    else
			lbms(i) = .false.
			hasbit = .true.
		    endif
490             continue
c               set hour to 00Z
                kpds(11) = 0
c		set unit to hours
                kpds(13) = 1

	        if (kpds(16).eq.3) then
c		    average values
		    kpds(14) = 0
		    kpds(15) = 24
		    kpds(17) = nsum
c                   kpds(20) = 4 - nsum
		else if (kpds(16).eq.10) then
c		    forecast or analysis
	            kpds(16) = 113
		    kpds(15) = 6
		    kpds(17) = nsum
c                   kpds(20) = 4 - nsum
		else
		    write(*,*) 'unexpected Time Range in pds'
		    write(*,*) 'need to extend program'
		    call exit(8)
		endif

		call setbit(kpds,kgds,hasbit)

                call ezgbwk(mean,lbms,kpds,kgds,1)
                icount = icount + 1
            else
                write(*,*) 'did not read - not enough data:',time2
            endif
        goto 250

500	continue

        write(*,*) 'no. of daily fields writen ',icount

        if (icount.eq.0) call exit(8)
        call exit(0)
        stop
        end