! koeppen.f90
! module procedure written by Tim Mitchell on 02.10.03

module Koeppen

implicit none

contains

!*******************************************************************************

subroutine GetKoeppen (Tmp,Pre,Boreal,KopNumb,KopAcro,KopName)

logical,dimension (12)		:: Winter,Summer
real,dimension (12),intent(in)	:: Tmp,Pre	! monthly tmp(degC),pre(mm)

logical,intent(in)		:: Boreal	! .TRUE. if boreal, else austral
integer,intent(out)		:: KopNumb	! Koeppen class as a number
character (len=3),intent(out)	:: KopAcro	! Koeppen class as an acronym
character (len=30),intent(out)	:: KopName	! Koeppen class as a name

real, parameter 	:: MissVal = -999.0

real :: TmpMin,PreMin,TmpMax,PreMax
real :: WinterPre,SummerPre,WinterPreMin,SummerPreMin,WinterPreMax,SummerPreMax
real :: Desert,Steppe,AnnualPre,AnnualTmp

integer :: Kop,XMonth

Winter=.FALSE. ; Summer=.FALSE.

if (Boreal) then
  Summer(4)=.TRUE. ; Summer(5)=.TRUE. ; Summer(6)=.TRUE. 
  Summer(7)=.TRUE. ; Summer(8)=.TRUE. ; Summer(9)=.TRUE. 
  Winter(1)=.TRUE. ; Winter(2)=.TRUE. ; Winter(3)=.TRUE. 
  Winter(10)=.TRUE. ; Winter(11)=.TRUE. ; Winter(12)=.TRUE. 
else
  Winter(4)=.TRUE. ; Winter(5)=.TRUE. ; Winter(6)=.TRUE. 
  Winter(7)=.TRUE. ; Winter(8)=.TRUE. ; Winter(9)=.TRUE. 
  Summer(1)=.TRUE. ; Summer(2)=.TRUE. ; Summer(3)=.TRUE. 
  Summer(10)=.TRUE. ; Summer(11)=.TRUE. ; Summer(12)=.TRUE. 
end if

! TmpMin=minval(Tmp) ; PreMin=minval(Pre) ; TmpMax=maxval(Tmp) ; PreMax=maxval(Pre)
! WinterPre=sum(Pre,Winter) ; SummerPre=sum(Pre,Summer)
! WinterPreMin=minval(Pre,Winter) ; SummerPreMin=minval(Pre,Summer)
! WinterPreMax=maxval(Pre,Winter) ; SummerPreMax=maxval(Pre,Summer)
! AnnualPre=sum(Pre) ; AnnualTmp=sum(Tmp)/12.0

TmpMin=100000 ; TmpMax=-100000 ; PreMin=100000 ; PreMax=-100000
AnnualTmp=0 ; AnnualPre=0 ; WinterPre=0 ; SummerPre=0
WinterPreMin=100000 ; WinterPreMax=-100000 ; SummerPreMin=100000 ; SummerPreMax=-100000
do XMonth=1,12
  if (Tmp(XMonth).LT.TmpMin) TmpMin=Tmp(XMonth)
  if (Tmp(XMonth).GT.TmpMax) TmpMax=Tmp(XMonth)
  if (Pre(XMonth).LT.PreMin) PreMin=Pre(XMonth)
  if (Pre(XMonth).GT.PreMax) PreMax=Pre(XMonth)
  
  AnnualTmp=AnnualTmp+Tmp(XMonth)
  AnnualPre=AnnualPre+Pre(XMonth)
  
  if (Winter(XMonth)) then
    WinterPre=WinterPre+Pre(XMonth)
    if (Pre(XMonth).LT.WinterPreMin) WinterPreMin=Pre(XMonth)
    if (Pre(XMonth).GT.WinterPreMax) WinterPreMax=Pre(XMonth)
  else
    SummerPre=SummerPre+Pre(XMonth)
    if (Pre(XMonth).LT.SummerPreMin) SummerPreMin=Pre(XMonth)
    if (Pre(XMonth).GT.SummerPreMax) SummerPreMax=Pre(XMonth)
  end if
end do
AnnualTmp=AnnualTmp/12.0

if (TmpMin.EQ.MissVal.OR.PreMin.EQ.MissVal) then	! data missing
  Kop=MissVal ; KopAcro="" ; KopName="unclassifiable"
else
  if      (WinterPre.GE.(0.7*AnnualPre)) then
    Steppe = AnnualTmp * 20.0
    Desert = AnnualTmp * 10.0
  else if (SummerPre.GE.(0.7*AnnualPre)) then
    Steppe = (AnnualTmp * 20.0) + 280
    Desert = (AnnualTmp * 10.0) + 140
  else
    Steppe = (AnnualTmp * 20.0) + 140
    Desert = (AnnualTmp * 10.0) +  70
  end if
  
  if      (AnnualPre.LE.Steppe) then			! Dry
    if (AnnualPre.LE.Desert) then				! Dry desert
      Kop=1 ; KopAcro="Bw" ; KopName="desert"
    else							! Dry steppe
      Kop=9 ; KopAcro="Bs" ; KopName="steppe"
    end if
    
    if (AnnualTmp.GT.18.0) then					! Dry warm 
      Kop=Kop+0	; KopAcro(2:2)="h" ; KopName="Dry warm "//KopName(1:6)		
    else							! Dry cool 
      Kop=Kop+4	; KopAcro(2:2)="k" ; KopName="Dry cool "//KopName(1:6)			
    end if
  else if (TmpMin.GT.18.0) then				! Tropical
    if      (PreMin.GT.60.0) then
      Kop=17 ; KopAcro="Af" ; KopName="Tropical always-wet"	! Tropical always-wet
    else if (AnnualPre.LT.(2540.0-(25.0*PreMin))) then
      Kop=21 ; KopAcro="Aw" ; KopName="Tropical wet-and-dry"	! Tropical wet-and-dry
    else							
      Kop=25 ; KopAcro="Am" ; KopName="Tropical monsoon"	! Tropical monsoon
    end if
  else if (TmpMax.GT.10.0) then	
    if (TmpMin.GT.-3.0) then
      Kop=29 ; KopAcro="C" ; KopName="Temp"		! Temperate
    else
      Kop=61 ; KopAcro="D" ; KopName="Cold"		! Cold
    end if
    
    if      (SummerPreMin.LT.40.0.AND.SummerPreMin.LT.(0.333*WinterPreMax)) then
      Kop=Kop+0	 ; KopAcro(2:2)="s" ; KopName=KopName(1:4)//" dry-summer"	! summer-dry 
    else if (WinterPreMin.LT.(0.1*SummerPreMax)) then
      Kop=Kop+12 ; KopAcro(2:2)="w" ; KopName=KopName(1:4)//" dry-winter"	! winter-dry 
    else
      Kop=Kop+24 ; KopAcro(2:2)="f" ; KopName=KopName(1:4)//" always-wet"	! always-wet 
    end if
    
    if      (TmpMax.GT.22.0) then
      Kop=Kop+0	; KopAcro(3:3)="a"  				! Temp/Cold warm-summer
      if (KopAcro(2:2).EQ."s") then
        KopName=KopName(1:4)//" warm-dry-summer"
      else
        KopName=KopName(1:15)//" warm-summer"
      end if
    else if (TmpMin.GT.-3.0.OR.count(Tmp.GT.10.0).GE.4) then
      Kop=Kop+4 ; KopAcro(3:3)="b"  				! Temp/Cold cool-summer
      if (KopAcro(2:2).EQ."s") then
        KopName=KopName(1:4)//" cool-dry-summer"
      else
        KopName=KopName(1:15)//" cool-summer"
      end if
    else
      Kop=Kop+8	; KopAcro(3:3)="c"  				! Temp/Cold brief-ssummer
      if (KopAcro(2:2).EQ."s") then
        KopName=KopName(1:4)//" brief-dry-summer"
      else
        KopName=KopName(1:15)//" brief-summer"
      end if
    end if
  else							! Polar
    if (TmpMax.GT.0) then
      Kop=97 ; KopAcro="Et" ; KopName="Polar tundra"		! Polar tundra
    else
      Kop=99 ; KopAcro="Ef" ; KopName="Polar permafrost"	! Polar permafrost
    end if
  end if
end if

KopNumb=MissVal
if (Kop.NE.MissVal) KopNumb=100-Kop

end subroutine GetKoeppen

!*******************************************************************************
! ths was originally used
! it was then modified into the form in GetKoeppen when it was found that on the
!   0.5deg grid some categories were only used for very small percentage of the
!   land surface; to simplify the classification these were removed

subroutine GetKoeppenOld (Tmp,Pre,Boreal,KopNumb,KopAcro,KopName)

logical,dimension (12)		:: Winter,Summer
real,dimension (12),intent(in)	:: Tmp,Pre	! monthly tmp(degC),pre(mm)

logical,intent(in)		:: Boreal	! .TRUE. if boreal, else austral
integer,intent(out)		:: KopNumb	! Koeppen class as a number
character (len=3),intent(out)	:: KopAcro	! Koeppen class as an acronym
character (len=30),intent(out)	:: KopName	! Koeppen class as a name

real, parameter 	:: MissVal = -999.0

real :: TmpMin,PreMin,TmpMax,PreMax
real :: WinterPre,SummerPre,WinterPreMin,SummerPreMin,WinterPreMax,SummerPreMax
real :: Desert,Steppe,AnnualPre,AnnualTmp

integer :: Kop,XMonth

Winter=.FALSE. ; Summer=.FALSE.

if (Boreal) then
  Summer(4)=.TRUE. ; Summer(5)=.TRUE. ; Summer(6)=.TRUE. 
  Summer(7)=.TRUE. ; Summer(8)=.TRUE. ; Summer(9)=.TRUE. 
  Winter(1)=.TRUE. ; Winter(2)=.TRUE. ; Winter(3)=.TRUE. 
  Winter(10)=.TRUE. ; Winter(11)=.TRUE. ; Winter(12)=.TRUE. 
else
  Winter(4)=.TRUE. ; Winter(5)=.TRUE. ; Winter(6)=.TRUE. 
  Winter(7)=.TRUE. ; Winter(8)=.TRUE. ; Winter(9)=.TRUE. 
  Summer(1)=.TRUE. ; Summer(2)=.TRUE. ; Summer(3)=.TRUE. 
  Summer(10)=.TRUE. ; Summer(11)=.TRUE. ; Summer(12)=.TRUE. 
end if

! TmpMin=minval(Tmp) ; PreMin=minval(Pre) ; TmpMax=maxval(Tmp) ; PreMax=maxval(Pre)
! WinterPre=sum(Pre,Winter) ; SummerPre=sum(Pre,Summer)
! WinterPreMin=minval(Pre,Winter) ; SummerPreMin=minval(Pre,Summer)
! WinterPreMax=maxval(Pre,Winter) ; SummerPreMax=maxval(Pre,Summer)
! AnnualPre=sum(Pre) ; AnnualTmp=sum(Tmp)/12.0

TmpMin=100000 ; TmpMax=-100000 ; PreMin=100000 ; PreMax=-100000
AnnualTmp=0 ; AnnualPre=0 ; WinterPre=0 ; SummerPre=0
WinterPreMin=100000 ; WinterPreMax=-100000 ; SummerPreMin=100000 ; SummerPreMax=-100000
do XMonth=1,12
  if (Tmp(XMonth).LT.TmpMin) TmpMin=Tmp(XMonth)
  if (Tmp(XMonth).GT.TmpMax) TmpMax=Tmp(XMonth)
  if (Pre(XMonth).LT.PreMin) PreMin=Pre(XMonth)
  if (Pre(XMonth).GT.PreMax) PreMax=Pre(XMonth)
  
  AnnualTmp=AnnualTmp+Tmp(XMonth)
  AnnualPre=AnnualPre+Pre(XMonth)
  
  if (Winter(XMonth)) then
    WinterPre=WinterPre+Pre(XMonth)
    if (Pre(XMonth).LT.WinterPreMin) WinterPreMin=Pre(XMonth)
    if (Pre(XMonth).GT.WinterPreMax) WinterPreMax=Pre(XMonth)
  else
    SummerPre=SummerPre+Pre(XMonth)
    if (Pre(XMonth).LT.SummerPreMin) SummerPreMin=Pre(XMonth)
    if (Pre(XMonth).GT.SummerPreMax) SummerPreMax=Pre(XMonth)
  end if
end do
AnnualTmp=AnnualTmp/12.0

if (TmpMin.EQ.MissVal.OR.PreMin.EQ.MissVal) then	! data missing
  Kop=MissVal ; KopAcro="" ; KopName="unclassifiable"
else
  if      (WinterPre.GE.(0.7*AnnualPre)) then
    Steppe = AnnualTmp * 20.0
    Desert = AnnualTmp * 10.0
  else if (SummerPre.GE.(0.7*AnnualPre)) then
    Steppe = (AnnualTmp * 20.0) + 280
    Desert = (AnnualTmp * 10.0) + 140
  else
    Steppe = (AnnualTmp * 20.0) + 140
    Desert = (AnnualTmp * 10.0) +  70
  end if
  
  if      (AnnualPre.LE.Steppe) then			! Dry
    if (AnnualPre.LE.Desert) then				! Dry desert
      Kop=1 ; KopAcro="Bw" ; KopName="desert"
    else							! Dry steppe
      Kop=9 ; KopAcro="Bs" ; KopName="steppe"
    end if
    
    if (AnnualTmp.GT.18.0) then					! Dry warm 
      Kop=Kop+0	; KopAcro(2:2)="h" ; KopName="Dry warm "//KopName(1:6)		
    else							! Dry cool 
      Kop=Kop+4	; KopAcro(2:2)="k" ; KopName="Dry cool "//KopName(1:6)			
    end if
  else if (TmpMin.GT.18.0) then				! Tropical
    if      (PreMin.GT.60.0) then
      Kop=20 ; KopAcro="Af" ; KopName="Tropical always-wet"	! Tropical always-wet
    else if (AnnualPre.LT.(2540.0-(25.0*PreMin))) then
      Kop=25 ; KopAcro="Aw" ; KopName="Tropical wet-and-dry"	! Tropical wet-and-dry
    else							
      Kop=30 ; KopAcro="Am" ; KopName="Tropical monsoon"	! Tropical monsoon
    end if
  else if (TmpMax.GT.10.0) then	
    if (TmpMin.GT.-3.0) then
      Kop=36 ; KopAcro="C" ; KopName="Temp"		! Temperate
    else
      Kop=60 ; KopAcro="D" ; KopName="Cold"		! Cold
    end if
    
    if      (SummerPreMin.LT.40.0.AND.SummerPreMin.LT.(0.333*WinterPreMax)) then
      Kop=Kop+0	 ; KopAcro(2:2)="s" ; KopName=KopName(1:4)//" dry-summer"	! summer-dry 
    else if (WinterPreMin.LT.(0.1*SummerPreMax)) then
      Kop=Kop+8	 ; KopAcro(2:2)="w" ; KopName=KopName(1:4)//" dry-winter"	! winter-dry 
    else
      Kop=Kop+16 ; KopAcro(2:2)="f" ; KopName=KopName(1:4)//" always-wet"	! always-wet 
    end if
    
    if      (TmpMax.GT.22.0) then
      Kop=Kop+0	; KopAcro(3:3)="a"  				! Temp/Cold warm-summer
      if (KopAcro(2:2).EQ."s") then
        KopName=KopName(1:4)//" warm-dry-summer"
      else
        KopName=KopName(1:15)//" warm-summer"
      end if
    else if (count(Tmp.GT.10.0).GE.4) then
      Kop=Kop+2 ; KopAcro(3:3)="b"  				! Temp/Cold cool-summer
      if (KopAcro(2:2).EQ."s") then
        KopName=KopName(1:4)//" cool-dry-summer"
      else
        KopName=KopName(1:15)//" cool-summer"
      end if
    else if (TmpMin.GT.-38.0) then
      Kop=Kop+4	; KopAcro(3:3)="c"  				! Temp/Cold brief-ssummer
      if (KopAcro(2:2).EQ."s") then
        KopName=KopName(1:4)//" brief-dry-summer"
      else
        KopName=KopName(1:15)//" brief-summer"
      end if
    else
      Kop=Kop+6	; KopAcro(3:3)="d"  				! Cold 	    hard-winter
      if (KopAcro(2:2).EQ."w") then
        KopName=KopName(1:4)//" hard-dry-winter"
      else
        KopName=KopName(1:15)//" hard-winter"
      end if
    end if
  else							! Polar
    if (TmpMax.GT.0) then
      Kop=93 ; KopAcro="Et" ; KopName="Polar tundra"		! Polar tundra
    else
      Kop=97 ; KopAcro="Ef" ; KopName="Polar permafrost"	! Polar permafrost
    end if
    
    if (TmpMin.LT.-38.0) then					! Polar hard-winter
      Kop=Kop+2	; KopAcro(3:3)="h" ; KopName=trim(KopName)//" hard-winter"			
    end if
  end if
end if

KopNumb=MissVal
if (Kop.NE.MissVal) KopNumb=100-Kop

end subroutine GetKoeppenOld

!*******************************************************************************

end module Koeppen