C........................................................Dec. 2021
C 1765 Kp storms 1994-2022 exclude kp peak at t0-24 : t0 at Kp peak 
C at least 2 days between the preceeding storm onset
C........................................................Apr.2021
C Extract Kp-index for -24h<t0< +48h for Kpmax >=5.0 events 
C 
C Select period of 72hrs starting from t0-24h, ending t0+48h; t0 at Kp peak 
C 
      PROGRAM decod5kpd1
C Input HSSUT.txt 1996-2010...........................Jan 2020
C add: select from -24h < to <+48h for HSS-events
C
C----------------------------------------------------------------
C
	DIMENSION inds(0:71) ! 3h indices for 24h (0r 72h file for 1h indices)
	+,ind24(0:23),ind8(0:7)
C
	CHARACTER*4 AYEAR
      CHARACTER*128 INFILE,INFHSS
	+,outfilkp1,outstrmkp1,outfilkp2,outstrmkp2
C     +,outfildst  !
	CHARACTER(10) DD
      CHARACTER(10) TT
      CHARACTER(5) ZZ
	CHARACTER*2 AYR,AMN,ADY
	CHARACTER*6 YYMMDD1,YYMMDD3
C++
C
      ICALLS=0
	CALL DATE_AND_TIME(DATE=DD,TIME=TT,ZONE=ZZ)
      TT(5:)='      '
      WRITE(*,*) 'PC Date: Year,Month,Day = ',DD,'Time = ',TT
C
	infile='d:/apkp/kpap/kpYEAR' !Kp
C
	outfilkp1='kpstrm3d1.txt' !total file for 1994-2022 superposed epoch moderate
	outfilkp2='kpstrm3d2.txt' !total file for 1994-2022 superposed epoch intense
	outstrmkp1='kpstorm1.txt'  ! moderate
	outstrmkp2='kpstorm2.txt'  ! intense
C
	nout=0
	OPEN(11,FILE=OUTFILKP1) ! superposed kp profile
	OPEN(12,FILE=OUTFILKP2)	! superposed kp profile
	OPEN(13,FILE=OUTSTRMKP1)   ! Catalog
	OPEN(14,FILE=OUTSTRMKP2)   ! Catalog
C
C Kp <1765 storms from 1994 to 2022
C
   28 format(1X,F8.3,3(1X,I2),3(1X,F6.0))
C   29 format(9X,2I2,2(1X,I2),5X,2I2,2(1X,I2),1X,I4,6X,2I2,2(1X,I2))
   29 format(1X,I4,1X,A6,1X,I2,2X,3(I2),1X,I2,1X,F4.1,1X,A6,1X,I2)
   30	format(1X,I4,1X,A6,1X,I2,2X,3(A2),1X,I2,1X,F4.1,1X,A6,1X,I2)
C Catalogue
  	infhss='kpstorm0.txt'  ! Kp storms
	open(10,file=infhss,action='READ')
C	read(10,*) txt
C
C Cycle on events ====================================================
C
      lda_pre=10
	iyear_pre=1995
	nout1=0		 ! moderate
	nout2=0		 ! intense
C	DO 177 iev=1,490   !
	DO 177 iev=1,1765   !
	iflag=0
C      read(10,29,end=27) iyr1,imn1,idy1,iut1,iyr2,imn2,idy2,iut2,idst
      read(10,29,end=27) nkp,YYMMDD1,iut1,iyr2,imn2,idy2,iut2,xkp
	+,YYMMDD3,iut3
	 if (xkp.lt.5.0) GOTO 177  ! ========>>>>>>>>>> Kp< 5.0
C?	iutmax=iut2/3	! in 3h intervals
C
C Kp peak:
C>>>>>>>>>>>>>>>>>>>>>>        iut2=iut2+1  ! middle hr at Kpmax in 3h interval
	call blet2(iyr2,AYR)
	call blet2(imn2,AMN)
	call blet2(idy2,ADY)
	iyear=1900+iyr2
	if (iyr2.lt.30) iyear=2000+iyr2
	 ryear=float(iyear)
         IF (xkp.lt.7.0) THEN
C   Moderate storm
	iflag=1
C
	                ELSE
C   Intense storm
	iflag=2
C
                             ENDIF
     	iutcnt=0
	do k=0,71
	inds(k)=0
	enddo
C
C
       lda1=ndoy(iyr2,imn2,idy2) ! day of Kp peak
C
C	lda1=ndoy(iyr1,imn1,idy1) ! 1d
      lda0=lda1-1 ! -1d
	lda2=lda1+2
C
      if ((iyear_pre.eq.iyear).and.(lda1-lda_pre.lt.2)) then
	iyear_pre=iyear
	lda_pre=lda1
	GOTO 177
	endif
	iyear_pre=iyear
	lda_pre=lda1
	ryear1=float(iyear1)
		if(int(ryear1/4.)*4.eq.iyear1) THEN 
	ndnyr=366
	                          ELSE
	ndnyr=365
	ENDIF
C
   37 format(12X,8(I2))	 ! Kp
   38 format(12X,8(I2),I3)	 ! Kp,Kpsum
C
	if (lda0.eq.0) then  !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
C
C read last day of prec. year:
C
	iyear0=iyear-1
       call blet4(iyear0,AYEAR)
	infile(16:19)=AYEAR
	open(15,file=infile,action='READ')
	ryear0=float(iyear0)
		if(int(ryear0/4.)*4.eq.iyear0) THEN 
	ndnyr=366
	                          ELSE
	ndnyr=365
	ENDIF
C
	DO 101 kk=1,ndnyr
	read(15,37,end=102) (ind8(ii),ii=0,7)
	ki=0
	do 99 k=0,7
C?	if (k.lt.iutmax) goto 99
	ind24(ki)=ind8(k)
	ind24(ki+1)=ind8(k)
	ind24(ki+2)=ind8(k)
	ki=ki+3
   99	continue
  101	CONTINUE
  102	close(unit=15)
C last day of prec year from iut1 to 23 h:   
       DO ii=iut2,23			! Corr!!!
       inds(iutcnt)=ind24(ii)
      iutcnt=iutcnt+1
	ENDDO
	endif				 !>>>>>>>>>>>>>>>>>>>>>>>>
C
C Prec. day from current year:
		if(int(ryear/4.)*4.eq.iyear) THEN 
	ndnyr=366
	                          ELSE
	ndnyr=365
	ENDIF
       call blet4(iyear,AYEAR)
	infile(16:19)=AYEAR
	open(15,file=infile,action='READ')
C
	if (lda0.gt.0) then !>>>>>>>>>>>>>>>>>>>>>>>>>>>
	DO 109 kk=1,lda0
	read(15,37,end=103) (ind8(ii),ii=0,7)
	ki=0
	do 98 k=0,7
C......	if (k.lt.iutmax) goto 99
	ind24(ki)=ind8(k)
	ind24(ki+1)=ind8(k)
	ind24(ki+2)=ind8(k)
	ki=ki+3
   98	continue
  109	CONTINUE
       DO ii=iut2,23		 ! Corr!!!
       inds(iutcnt)=ind24(ii)
      iutcnt=iutcnt+1
	ENDDO
	endif				 !>>>>>>>>>>>>>>>>>>>>>>>>
C 1st day from current year:
	read(15,38,end=103) (ind8(ii),ii=0,7),kpsum
	ki=0
	do k=0,7
	ind24(ki)=ind8(k)
	ind24(ki+1)=ind8(k)
	ind24(ki+2)=ind8(k)
	ki=ki+3
	enddo
C	  1st day :
       DO ii=0,23
       inds(iutcnt)=ind24(ii)
      iutcnt=iutcnt+1
	ENDDO
C 2nd day from current year:
	read(15,37,end=103) (ind8(ii),ii=0,7)
	ki=0
	do k=0,7
	ind24(ki)=ind8(k)
	ind24(ki+1)=ind8(k)
	ind24(ki+2)=ind8(k)
	ki=ki+3
	enddo
C	  2nd day :
       DO ii=0,23
       inds(iutcnt)=ind24(ii)
      iutcnt=iutcnt+1
	ENDDO
C 3rd day from current year:
	read(15,37,end=103) (ind8(ii),ii=0,7)
	ki=0
	do k=0,7
	ind24(ki)=ind8(k)
	ind24(ki+1)=ind8(k)
	ind24(ki+2)=ind8(k)
	ki=ki+3
	enddo
C	  3rd day :
       DO ii=0,23
	if (iutcnt.lt.72) then
       inds(iutcnt)=ind24(ii)
      iutcnt=iutcnt+1
              	else
	exit
              	endif
	ENDDO
  103	close(unit=15)
C
C Next year if necessary:
       if (iutcnt.lt.72) then   !>>>>>>>>>>>>>>>>>>>>>>>>>>>
	iyear1=iyear+1
	call blet4(iyear1,AYEAR)
	infile(16:19)=AYEAR
	open(15,file=infile,action='READ')
C
	read(15,37,end=104) (ind8(ii),ii=0,7)
	ki=0
	do k=0,7
	ind24(ki)=ind8(k)
	ind24(ki+1)=ind8(k)
	ind24(ki+2)=ind8(k)
	ki=ki+3
	enddo
C	  next day :
       DO ii=0,23
	if (iutcnt.lt.72) then
       inds(iutcnt)=ind24(ii)
      iutcnt=iutcnt+1
              	else
	exit
              	endif
	ENDDO
  104	close(unit=15)
      endif	                !>>>>>>>>>>>>>>>>>>>>>
C
  138 format(72(1X,I2))
  139	FORMAT(1X,F8.3,1X,I3,2(1X,A2),1X,I2,1X,F3.1,2(1X,I3))
C Check Kp max at t0-24 : t0 Exclude these results
      maxkp=0
	do k=0,71
	if (inds(k).gt.maxkp) then
	k1=k
	maxkp=inds(k)
	endif
	enddo
	ykpmax=maxkp/10.
Crem	if (k1.le.23) goto 177  ! avoid output for Kpmax at t0-24h : t0
	if ((k1.le.23).or.(k1.gt.47).or.(ykpmax.gt.xkp)) goto 177
      xyear=ryear+(lda1-1.)/ndnyr
	iut2=iut2+1
C Output matrices:
        IF (iflag. eq.1) THEN
C output new Catalogue 5.0 <= Kp < 7.0
 	nout1=nout1+1
C
      write(*,139) xyear,lda1,amn,ady,iut2,xkp,kpsum,nout1
	write(13,139) xyear,lda1,amn,ady,iut2,xkp,kpsum,nout1
C
 	write(*,138) (inds(k),k=0,71)
	write(11,138) (inds(k),k=0,71)
	                 ELSE
C output new Catalogue Kp >=7.0
 	nout2=nout2+1
C
      write(*,139) xyear,lda1,amn,ady,iut2,xkp,kpsum,nout2
	write(14,139) xyear,lda1,amn,ady,iut2,xkp,kpsum,nout2
C
 	write(*,138) (inds(k),k=0,71)
	write(12,138) (inds(k),k=0,71)
	                     ENDIF
  177	CONTINUE
   27 close(unit=10)
      close(unit=11)
C
C
   35      write(*,*) infile
      pause ' ' 
      STOP
       END

C***********************************************************
