      program extwcpia2
C..Produce 1h North:dWn(pia) South: dWs(pia) indices..for 24 Mlat=+-60o sites..........July 2022
C Output Range; dW=Wcmax-WCmin
C JPL Wc indices 
C Wc interp from 4 W indices: at dlati=2.5o, dlong=5o around Kp station
C Wc = combined  with weights from 4 W-indices surrounding Kp station
C 
C 
C
      CHARACTER*64 outfilen,outfiles
	CHARACTER*2 AYR,AMN,ADY,AUT,aa(24),bb
	CHARACTER*4 AYEAR,AMEAN
	DIMENSION wn(24,0:23),ws(24,0:23),glatn(24),glonn(24)
	+,glats(24),glons(24),IM(12),iwn(0:23),iws(0:23)
	DATA im/31,28,31,30,31,30,31,31,30,31,30,31/
	DATA glatn/ 50.8,  51.0,  51.8,  53.0,  54.6,  56.6,  58.9,  61.3
	+,  63.7,  66.0,  67.8,  69.1,  69.5,  69.1,  67.8,  66.0,  63.7
     +,  59.6,  58.9,  56.6,  54.6,  53.0,  51.8,  51.0/
	DATA glats/-69.5, -69.1, -67.8, -66.0, -63.7, -61.3, -58.9, -56.6
	+, -54.6, -53.0, -51.8, -51.0, -50.8, -51.0, -51.8, -53.0, -54.6
     +, -56.6, -58.9, -61.3, -63.7, -66.0, -67.8, -69.1/
	DATA glonn/287.3, 299.1, 311.0, 323.1, 335.4, 348.2, 358.5,  15.7
	+,  31.0,  47.7,  66.1,  86.2, 107.3, 128.4, 148.5, 166.9, 176.3
     +, 199.0, 213.1, 226.5, 239.2, 251.5, 263.6, 275.5/
	DATA glons/287.3, 308.4, 328.5, 347.0,   3.7,  19.0,  33.1,  46.5
	+,  59.2,  71.5,  83.6,  94.5, 107.3, 119.1, 131.0, 143.1, 155.4
     +, 168.2, 181.5, 195.7, 211.0, 227.7, 246.1, 266.2/
	DATA aa/'00','01','02','03','04','05','06','07','08','09','10'
	+,'11','12','13','14','15','16','17','18','19','20','21','22','23'/
C
C
	bb='UT'
	amean='mean'
	outfilen='d:\web\graf\datw\wpoyrn.txt'	 ! daily 0:23 UT + daily Wpo index
	outfiles='d:\web\graf\datw\wpoyrs.txt'
C
198	format(4X,A2,24(4X,A2),3X,A5)
C
C Year cycle
       DO 150 iyear=1994,2022
	 call blet4(iyear,AYEAR)
	 AYR=Ayear(3:4)
	 outfilen(21:22)=AYR          
       OPEN(12,file=outfilen)
	 outfiles(21:22)=AYR          
       OPEN(13,file=outfiles)
	write(*,198) bb,aa,amean
	write(12,198) bb,aa,amean
      write(13,198) bb,aa,amean
C
	ndy=0
C
C++ READ 1h Wc
		z1=iyear/4.0
      jz=int(z1)*4
      IF(jz.EQ.iyear) THEN
	ndnr=366
	im(2)=29
		  ELSE
                ndnr=365
	          im(2)=28
		       ENDIF
C
       imn2=12
	if (iyear.eq.2022) imn=06  !TEMPPPPPPPPPPPPPPPP
C
C Month cycle
	 DO 151 imn=1,imn2
	call blet2(imn,amn)
	ndmn=im(imn)
C
C
C daily cycle
       DO 152 idy=1,ndmn 
	call blet2(idy,ady)
C
	lda=ndoy(iyr,imn,idy)
C
C UT cycle
  43    DO 153 iut=0,23
	call blet2(iut,aut)
	DO k=1,24
	wn(k,iut)=0.
	ws(k,iut)=0.
	ENDDO
C
C cycle on stations
C
       DO ista=1,24
C
      glat=glatn(ista)
	glon=glonn(ista)
	call subinterpw(glat,glon,ayr,amn,ady,aut,res)
	wn(ista,iut)=res
      glat=glats(ista)
	glon=glons(ista)
	call subinterpw(glat,glon,ayr,amn,ady,aut,res)
	ws(ista,iut)=res
	ENDDO
  153 CONTINUE  ! End UT cycle
C
        wnday=0.
	wsday=0.
         DO kut=0,23
	 wnmax=0.
	wnmin=0.
	 wsmax=0.
	wsmin=0.
C find max, min:
        do kk=1,24
	if (wnmax.lt.wn(kk,kut)) wnmax=	wn(kk,kut)
	if (wnmin.gt.wn(kk,kut)) wnmin=	wn(kk,kut)
	if (wsmax.lt.ws(kk,kut)) wsmax=	ws(kk,kut)
	if (wsmin.gt.ws(kk,kut)) wsmin=	wn(kk,kut)
	  enddo
	dwn=wnmax-wnmin
	dws=wsmax-wsmin
	iwn(kut)=nint(dwn*100.)
	iws(kut)=nint(dws*100.)
	wnday=wnday+dwn
	wsday=wsday+dws
	  ENDDO
	iwnday=nint(wnday/24.*100.)
	iwsday=nint(wsday/24.*100.)
C Output
      write(*,32) ayr,amn,ady,(iwn(k),k=0,23),iwnday
      write(12,32) ayr,amn,ady,(iwn(k),k=0,23),iwnday
      write(*,32) ayr,amn,ady,(iws(k),k=0,23),iwsday
      write(13,32) ayr,amn,ady,(iws(k),k=0,23),iwsday
  152 CONTINUE  ! End Day cycle
C
  151 CONTINUE  ! End Mon cycle
   50		close(unit=12)
      	close(unit=13)
  150 continue  ! End Year cycle
   32 format(3(A2),1X,24(1X,I5),2X,I5)
	pause ' '
	STOP
	END