#include "dims.h"
      SUBROUTINE CMPNOZ
      use cons_module,only: len1,len2,kmax,kmaxp1,rmassinv_o2
      implicit none
C     ****
C     ****     ADVANCES PSI(NOZ) BY ONE TIME STEP
C     ****
#include "params.h"
#include "fcom.h"
#include "vscr.h"
#include "index.h"
#include "buff.h"
#include "crates_const.h"
#include "crates_tdep.h"
#include "phys.h"
#include "lowbnd.h"
#include "compcom.h"
#include "mwt.h"
#include "diffk.h"
#include "cmpdat.h"
!
! Local:
      character(len=80) :: title
      real :: PHINO(3)=(/0.814,0.866,0.926/)
      real :: alfa,xynoz
      integer :: i,k,ibnd,ibndb
      integer :: npnozk,npn4sk,npnok,npno2k,npohk,npn2dk,ntk,nps1k,
     |  nps2k,npo1k,npo3k,npho2k,nmsk,nno2k,nek,no2pk,nphoxk
!
      NPNOZK = NJ+NPNOZ
      NPN4SK = NJ+NPN4S
      NPNOK = NJ+NPNO
      NPNO2K = NJ+NPNO2
      NPOHK = NJ+NPOH
      NPN2DK = NPN2D+NJNP
      NTK = NJ+NT
      NPS1K = NJ+NPS
      NPS2K = NPS1K+KMAXP1
      NPO1K = NJ+NPO1
      NPO3K = NJ+NPO3
      NPHO2K = NJ+NPHO2
      NMSK = NJ+NMS
      NNO2K = NNO2
      NEK = NJ+NE
      NO2PK = NJ+NO2P
      NPHOXK = NJ+NPHOX
C     ****
      DO 3 I=1,LEN2
C     ****
C     ****     S3 = N*MBAR   (K+1/2)
C     ****
      S3(I,1) = XNMBAR(I,1)
C     ****
C     ****     S4(I,1) = RMNOZ(TRUE)
C     ****
      S1(I,1) = (BETA9(I,1)*F(I,NPO3K)/RMO3+BETA10(I,1)*
     1          F(I,NPHOXK)/RMTRU(I,1)*RATIO2(I,1)*RATIO3(I,1))
     2          /(BETA11*F(I,NPO1K)/RMO1
     3          +BETA12(I,1)*F(I,NPO3K)/RMO3+0.5*(
     4          XJNO2(I,1)+XJNO2(I+LEN1,1))/S3(I,1)+BETA13*
     5          F(I,NPN4SK)/RMN4S)
      S2(I,1) = 1./(1.+S1(I,1))
      if (s2(i,1) < 1.e-6) s2(i,1) = 1.e-6
      if (s1(i,1) < (1.-s2(i,1))/s2(i,1)) s1(i,1) = (1.-s2(i,1))/s2(i,1)
      S4(I,1) = S2(I,1)*(RMNO+S1(I,1)*RMNO2)
      ratio4(i,1) = s1(i,1)
      if (s1(i,1) < 1.e-6) ratio4(i,1) = 1.e-6
      ratio5(i,1) = s2(i,1)
      if (s2(i,1) < 1.e-6) ratio5(i,1) = 1.e-6
      rmtru1(i,1) = s4(i,1)
      if (s4(i,1) < 1.e-6) rmtru1(i,1) = 1.e-6
C     RATIO4(I,1) = S1(I,1)
C     RATIO5(I,1) = S2(I,1)
C     RMTRU1(I,1) = S4(I,1)
  3   CONTINUE
C
      DO 4 I=1,LEN2

C     ****
C     ****     S2 = (NUMBER DENSITY PRODUCTION) * (RMNOZ(TRUE)/RMNO)
C     ****
      S2(I,1) = S3(I,1)**2*(F(I,NPS1K)*rmassinv_o2*(BETA1(I,1)*
     1          F(I,NPN4SK)/RMN4S+BETA2*F(I,NPN2DK)/
     2          RMN2D)+BETA8*F(I,NPHOXK)/RMTRU(I,1)
     3          *RATIO1(I,1)*RATIO3(I,1)*F(I,NPN4SK)/RMN4S)
     4          +2.*BETA14*XNO1D(I,1)/RMO1*S3(I,1)*
     5          RMN2O(I,1)*S3(I,1)/(.5*(F(I,NMSK)+F(I,NMSK+1)))
      S2(I,1) = S2(I,1)*RMTRU1(I,1)/RMNO
C     ****
C     ****     S1 = (NUMBER DENSITY LOSS)/N(NOZ)
C     ****
      S1(I,1) = -(S3(I,1)*(BETA3(I,1)*F(I,NPN4SK)/RMN4S+
     1          BETA6*F(I,NPN2DK)/RMN2D)+.5*(XJNO(I,1)+
     2          XJNO(I+LEN1,1)+XJNOP(I,1)+XJNOP(I+LEN1,1)+
     3          XJNOPN(I,1)+XJNOPN(I+LEN1,1)+DNOEUV(I,1)+
     4          DNOEUV(I+LEN1,1))+RK5
     5          *F(I,NO2PK))*RATIO5(I,1)-S3(I,1)*(BETA13*
     6          F(I,NPN4SK)/RMN4S+BETA12(I,1)*F(I,NPO3K)/RMO3)
     7          *RATIO4(I,1)*RATIO5(I,1)
  4   CONTINUE
c     write(title,"('RMNOZ PRODUCTION (S2 FROM CMPNOZ): LAT=',i2)") j
c     call mkcon(s2,zimxp,zkmxp,0.,0,'LONGITUDE','PRESSURE',title)
C     ****
      DO 5 I=1,LEN1
      T1(I) = 0.
      T2(I) = 1.
C     ****
C     ****     NUMBER DENSITY MIXING RATIO OF NOZ AT LOWER BOUNDARY
C     ****       RMOZ(TRUE) AT LOWER BOUNDARY
C     ****
      T3(I) = -XNOZLB(J)*RMTRU1(I,1)/F(I,NMSK)
C     ****     DIFFUSIVE EQUILIBRIUM AT UPPER BOUNDARY
      T4(I) = 0.
    5 CONTINUE
      IBND = 0
      IBNDB = 0
      ALFA=0.
      XYNOZ=1.E-3
      CALL MINOR(NPNOZ,NPNOZNM,RMNO,PHINO,ALFA,IBND,IBNDB,WNOZ,
     1XYNOZ,NPDHNOZ,difkk)
C     ****
C     **** PARTITION NOZ INTO NO2 AND NO
C     ****
      NPNOZK = NJNP+NPNOZ-1
      NPNOK = NJNP+NPNO-1
      NPNO2K = NJNP+NPNO2-1
      DO 6 K=1,KMAX
      NPNOZK = NPNOZK+1
      NPNOK = NPNOK+1
      NPNO2K = NPNO2K+1
      DO 6 I=1,LEN1
      S11(I,K) = XYNOZ*WNOZ(K)
      if (f(i,npnozk) < s11(i,k)) f(i,npnozk) = s11(i,k)
      F(I,NPNOK) = RATIO5(I,K)*F(I,NPNOZK)*RMNO/RMTRU1(I,K)
      F(I,NPNO2K) = RATIO4(I,K)*RATIO5(I,K)*F(I,NPNOZK)*RMNO2/
     1              RMTRU1(I,K)
  6   CONTINUE
C     ****
C     ****     SET LEVEL KMAX+1 USING LOGARITHMIC INTERPOLATION
C     ****
      NPNOK = NJNP+NPNO+KMAX
      NPNO2K = NJNP+NPNO2+KMAX
      DO I = 1,LEN1
        F(I,NPNOK) = F(I,NPNOK-1)**2/F(I,NPNOK-2)
        F(I,NPNO2K) = F(I,NPNO2K-1)**2/F(I,NPNO2K-2)
      ENDDO
      RETURN
      END
C