#include "dims.h"
      SUBROUTINE CPKTKM
      use cons_module,only: len1,len3,kmax,kmaxp1,t0,rmassinv,gask
      implicit none
C     ****     CALCULATES CP, KT AND KM
#include "params.h"
#include "fcom.h"
#include "vscr.h"
#include "index.h"
#include "buff.h"
#include "phys.h"
!
! Local:
      integer :: nps1k,nps2k,ncpk,nktk,nkmk,i,ntk,k
!

!     call addfsech('tn_cp',' ',' ',f(1,nj+nt)  ,zimxp,zkmxp,zkmxp,j)
!     call addfsech('o2_cp',' ',' ',f(1,nj+nps) ,zimxp,zkmxp,zkmxp,j)
!     call addfsech('o1_cp',' ',' ',f(1,nj+nps2),zimxp,zkmxp,zkmxp,j)

      NPS1K=NJ+NPS
      NPS2K=NPS1K+KMAXP1
      NCPK=NCP
      NKTK=NKT
      NKMK=NKM
      DO 1 I=1,LEN3
C     ****     S1=MBAR
      S1(I,1)=1./(F(I,NPS1K)*rmassinv(1)+F(I,NPS2K)*rmassinv(2)+
     1(1.-F(I,NPS1K)-F(I,NPS2K))*rmassinv(3))
C     ****     S2=PO2
      S2(I,1)=S1(I,1)*F(I,NPS1K)*rmassinv(1)
C     ****     S3=PO
      S3(I,1)=S1(I,1)*F(I,NPS2K)*rmassinv(2)
C     ****     S4=PN2
      S4(I,1)=1.-S2(I,1)-S3(I,1)
C     ****     CP=R/2.*(PO2*7./32.+PN2*7./28.+PO*5./16.)
      F(I,NCPK)=gask*.5*(S2(I,1)*7./32.+S4(I,1)*7./28.+S3(I,1)*5./16.)
C     ****     KM=(PO2*4.03+PN2*3.42+PO*3.9)
      F(I,NKMK)=S2(I,1)*4.03+S4(I,1)*3.42+S3(I,1)*3.9
C     ****     KT=(PO2+PN2)*56.+PO*75.9
      F(I,NKTK)=(S2(I,1)+S4(I,1))*56.+S3(I,1)*75.9
    1 CONTINUE

!     call addfsech('fcp1',' ',' ',f(1,ncp),zimxp,zkmxp,zkmxp,j)
!     call addfsech('fkm1',' ',' ',f(1,nkm),zimxp,zkmxp,zkmxp,j)
!     call addfsech('fkt1',' ',' ',f(1,nkt),zimxp,zkmxp,zkmxp,j)

      NTK=NJ+NT-1
      NKTK=NKTK-1
      NKMK=NKMK-1
      DO 2 K=1,KMAX
      NKTK=NKTK+1
      NKMK=NKMK+1
      NTK=NTK+1
      DO 2 I=1,LEN1
      F(I,NKMK)=F(I,NKMK)*(F(I,NTK)+.5*(T0(K)+T0(K+1)))**0.69*1.E-6
      F(I,NKTK)=F(I,NKTK)*(F(I,NTK)+.5*(T0(K)+T0(K+1)))**0.69
    2 CONTINUE
      NKMK=NKMK+1
      NKTK=NKTK+1
      DO 3 I=1,LEN1
      F(I,NKMK)=1.E-6*F(I,NKMK)*(F(I,NTK)+1.5*T0(KMAXP1)-.5*T0(KMAX))
     1**0.69
      F(I,NKTK)=F(I,NKTK)*(F(I,NTK)+1.5*T0(KMAXP1)-.5*T0(KMAX))**0.69
    3 CONTINUE

!     call addfsech('fcp2',' ',' ',f(1,ncp),zimxp,zkmxp,zkmxp,j)
!     call addfsech('fkm2',' ',' ',f(1,nkm),zimxp,zkmxp,zkmxp,j)
!     call addfsech('fkt2',' ',' ',f(1,nkt),zimxp,zkmxp,zkmxp,j)

      NKMK=NKMK+1
      NKTK=NKTK+1
      NCPK=NCPK+KMAXP1
      DO 4 K=1,KMAX
      NKMK=NKMK-1
      NKTK=NKTK-1
      NCPK=NCPK-1
      DO 4 I=1,LEN1
      F(I,NCPK)=.5*(F(I,NCPK)+F(I,NCPK-1))
      F(I,NKMK)=.5*(F(I,NKMK)+F(I,NKMK-1))
      F(I,NKTK)=.5*(F(I,NKTK)+F(I,NKTK-1))
    4 CONTINUE
      DO 5 I=1,LEN1
      F(I,NCPK-1)=2.*F(I,NCPK-1)-F(I,NCPK)
      F(I,NKMK-1)=2.*F(I,NKMK-1)-F(I,NKMK)
      F(I,NKTK-1)=2.*F(I,NKTK-1)-F(I,NKTK)
    5 CONTINUE

!     call addfsech('CP',' ',' ',f(1,ncp),zimxp,zkmxp,zkmxp,j)
!     call addfsech('KT',' ',' ',f(1,nkt),zimxp,zkmxp,zkmxp,j)
!     call addfsech('KM',' ',' ',f(1,nkm),zimxp,zkmxp,zkmxp,j)

      RETURN
      END
C