      program extlatWs
C.........................................................May 2021
C Shub W(GTEC) index
C.........................................................Dec. 2020
C Include weight e(glat)
C
C.........................................................Nov. 2020
C Weighted PIA -pos ion activity, NIA - neg ion activity
C extracted from global Wc(TEC) index derived from JPL GIM-TEC (IONEX)
C for ADY_pre and forecast for DD2
C
	DIMENSION weg(4),IM(12),iwneg(35),iwpos(35),iwran(35) ! LATS 85, 80, 75,...,0,...,85  :dfi=5deg
	CHARACTER*128 INFILE,OUTFILEN,OUTFILEP,OUTFILER !global weighted Wfp,Wfn,Wfr = range index
	CHARACTER*2 AYR,AUT,AMN,ADY,ADY1,CLAT
	CHARACTER*4 AYEAR
	CHARACTER*3 ALAT,AUTOUT(35)
	INTEGER ivsr1(0:72),ivsr2(0:72) !  set of Wf indices for long=-180:5:180 at given lati
		DATA WEG/0.,0.6931,1.0986,1.3863/  ! =ln1, ln2, ln3, ln4
            DATA IM/31,28,31,30,31,30,31,31,30,31,30,31/
C
C
      infile='d:\web\gtc\YY\MM\wgMMDDUT.sYY'   ! current mn
      outfilen='d:\web\graf\dats\wlnYYslat.txt'   !  Cnt(W-) 	 WL
      outfilep='d:\web\graf\dats\wupYYslat.txt'   !  Cnt(W+) 	 WU
      outfiler='d:\web\graf\dats\werYYslat.txt'   !  Cnt(Wrange) WE
C
	iglat=-90
      do ilat=1,35
	iglat=iglat+5
	iglat2=abs(iglat)
	call blet2(iglat2,clat)
	ALAT='   '
	ALAT(2:3)=clat
	if (iglat.lt.0) ALAT(1:1)='-'
	autout(ilat)=alat
	enddo
   17 format(A2)
   31 FORMAT('YYMMDD',1X,35(1X,A3,1X))
   32 FORMAT(73(1X,I4))
C
        write(*,*) ' Enter iyear:'
	read(*,*) iyear
        write(*,*) ' Enter imn:'
	read(*,*) imn
        write(*,*) ' Enter idy1,idy2:'
	read(*,*) idy1,idy2
	 
Crem         DO 700 iyear=1994,2020 ! YEARRRRRRRRRRRRRRRRRRRRRRRRRR
	call blet4(iyear,AYEAR)
	IYR=iyear-iyear/100*100   ! current yr
	z1=iyear/4.0
      jz=int(z1)*4

      IF(jz.EQ.iyear) THEN
               IM(2)=29
	dnr=366.
        ELSE
                IM(2)=28
	  	dnr=365.
	       ENDIF
      AYR=AYEAR(3:4)  ! for current day  	
	infile(12:13)=AYR
	infile(28:29)=AYR
	outfilen(21:22)=AYR
	outfilep(21:22)=AYR
	outfiler(21:22)=AYR
	OPEN(11,file=outfilen,access='APPEND')
	OPEN(12,file=outfilep,access='APPEND')
	OPEN(13,file=outfiler,access='APPEND')
C
C
Crem	    DO 51 imn=1,12        ! MONNNNNNNNNNNNNNNNN
	call blet2(imn,AMN)
	infile(15:16)=AMN
	infile(20:21)=AMN
CREM         idy1=1
CREM	idy2=im(imn)
CREM	mndy=im(imn)
       call blet2(idy1,ADY1)
C
C  Day-to-day cycle
          DO 52 idy=idy1,idy2 !DAYYYYYYYYYYYYYY
	call blet2(idy,ADY)
	infile(22:23)=ADY
C
  330	CONTINUE
C
C	IF ((AMN.eq.'01').and.(ADY.eq.ADY1)) THEN
	IF (ADY.eq.ADY1) THEN
	write(*,31) (autout(k),k=1,35)
	write(11,31) (autout(k),k=1,35)
	write(12,31) (autout(k),k=1,35)
	write(13,31) (autout(k),k=1,35)
	ENDIF
C  UT cycle
C
      DO k=1,35
	iwneg(k)=0
	iwpos(k)=0
	iwran(k)=0
	ENDDO
	 wdyp=0.
	wdyn=0.
	wdyr=0.
C
C-      glat=-90.
         DO 100 iout=1,35
	sumn=0.		 ! count of Wf >=1
	sump=0.		 ! count of Wf <=-1
	icntp=0
	swegp=0.
      icntn=0
	swegn=0.
C++	iout=0
	il2=iout*2
	il1=il2-1
C-	glat=glat+2.5				  ! ADC
	DO 53 iut=0,23  ! UTTTTTTTTTTTTTT
	call blet2(iut,AUT)
      infile(24:25)=aut
C
   	OPEN(10,file=infile,action='READ')
Crem	sumn=0.		 ! count of Wf >=1
Crem	sump=0.		 ! count of Wf <=-1
Crem	icntp=0
Crem	swegp=0.
Crem      icntn=0
Crem	swegn=0.
Crem	iout=0
C
C/ cycle for lats:
       glat=-90.0	                  ! AD
	  do ilat=1,il1
	glat=glat+2.5
	  read(10,32) (ivsr1(k),k=0,72) 
	enddo
	read(10,32) (ivsr2(k),k=0,72) 
	close(unit=10)
C
	DO 333 ilon=0,72
          if (ivsr1(ilon).gt.1) then
	elat=1.0+0.0001235*glat**2	  ! AD
        iad=ivsr1(ilon)
c	sump=sump+weg(iad)*ivsr(ilon)
	sump=sump+elat*weg(iad)*ivsr1(ilon)  ! AD
	icntp=icntp+1
	swegp=swegp+weg(iad)
Crem	goto 333
	endif
          if (ivsr2(ilon).gt.1) then
	elat=1.0+0.0001235*(glat+2.5)**2	  ! AD
        iad=ivsr2(ilon)
c	sump=sump+weg(iad)*ivsr(ilon)
	sump=sump+elat*weg(iad)*ivsr2(ilon)  ! AD
	icntp=icntp+1
	swegp=swegp+weg(iad)
Crem	goto 333
	endif
          if (ivsr1(ilon).lt.-1) then
	elat=1.0+0.0001235*glat**2	  ! AD
        iad=-ivsr1(ilon)
c	sumn=sumn+weg(iad)*ivsr(ilon)
	sumn=sumn+elat*weg(iad)*ivsr1(ilon)   ! AD
	icntn=icntn+1
	swegn=swegn+weg(iad)
Crem	goto 333
	endif
          if (ivsr2(ilon).lt.-1) then
	elat=1.0+0.0001235*(glat+2.5)**2	  ! AD
        iad=-ivsr2(ilon)
c	sumn=sumn+weg(iad)*ivsr(ilon)
	sumn=sumn+elat*weg(iad)*ivsr2(ilon)   ! AD
	icntn=icntn+1
	swegn=swegn+weg(iad)
Crem	goto 333
	endif
  333	CONTINUE  ! end long cycle
C-   54 CONTINUE  ! end lat cycle
   53	CONTINUE           ! end UT cycle
      xpos=0.
      if (icntp.gt.0.) then
C		xpos=sump/swegp		 ! /sum(weights)
		xpos=sump/icntp		 ! /cnt(Wf+)
	endif
      xneg=0.
      if (icntn.gt.0.) then
C			xneg=sumn/swegn	  ! /sum(weights)
			xneg=sumn/icntn	  ! /cnt(Wf-)
	endif
	xrange=xpos-xneg
	wdyp=wdyp+xpos
	wdyn=wdyn+xneg
      wdyr=wdyr+xrange
Crem      iwneg(iut)=nint(xneg*100.)
Crem	iwpos(iut)=nint(xpos*100.)
Crem	iwran(iut)=nint(xrange*100.)
      iwneg(iout)=nint(xneg*100.)
	iwpos(iout)=nint(xpos*100.)
	iwran(iout)=nint(xrange*100.)
  100	CONTINUE         ! lats-out cycle
c//	iout=iout+1
C-       iwpos(24)=nint(wdyp/24.*100.)
C-	iwneg(24)=nint(wdyn/24.*100.)
C-	iwran(24)=nint(wdyr/24.*100.)
	write(*,33) AYR,AMN,ADY,(iwpos(k),k=1,35)
	write(12,33) AYR,AMN,ADY,(iwpos(k),k=1,35)
	write(*,33) AYR,AMN,ADY,(iwneg(k),k=1,35)
	write(11,33) AYR,AMN,ADY,(iwneg(k),k=1,35)
	write(*,33) AYR,AMN,ADY,(iwran(k),k=1,35)
	write(13,33) AYR,AMN,ADY,(iwran(k),k=1,35)
   33 format(3A2,35(1X,I4))
   52 CONTINUE  ! daily cycle
C//       call subws3g(AYR,AMN,idy1,idy2)
   51 CONTINUE  ! Monthly cycle
      close(unit=11)
	close(unit=12)
	close(unit=13)
  700 CONTINUE              ! END YEAR CYCLE
      pause ' '
	STOP
	END
