!-------------------------------------------------------------------------------
! GEN_CO2 — gerador offline das funcoes de transmissao de CO2 (GFDL LW) no
! layout lido pelo CONRAD do etafcst (CO2.dat / co2.${LM}_${PT}mb_${ppm}ppm).
!
! Uso:  gen_co2.x <arquivo_deta> <PT_hPa> <ppm> <arquivo_saida>
!   ex: gen_co2.x deta_80_1mb 1.0 330 co2.80_1mb_330ppm
! Requer no diretorio corrente: tr49t85, tr49t67, tr67t85 (tabelas 109x109
! Schwarzkopf-Fels, 3 perfis de T cada, formatadas F20.14 — ex.: WRF/run).
!
! O RATIO passado ao kernel = ppm/330 (tabelas-base sao 330 ppmv).
! Registros de saida (ordem exata do CONRAD.f90):
!   2 x LP1 (STEMP,GTEMP); 6 x LM (CDTM51,CO2M51,C2DM51,CDTM58,CO2M58,C2DM58);
!   6 x LP1*LP1 (CDT51,CO251,C2D51,CDT58,CO258,C2D58) com I2 variando mais
!   rapido p/ I1 fixo; 6 x LP1 (CDT31,CO231,C2D31,CDT38,CO238,C2D38);
!   6 x LP1 (CDT71,CO271,C2D71,CDT78,CO278,C2D78).
!-------------------------------------------------------------------------------
      PROGRAM GEN_CO2
      USE GFDLCO2
      IMPLICIT NONE
      CHARACTER(256) :: detafile, outfile, arg
      CHARACTER(16)  :: mode
      REAL           :: pt_hpa, ppm, pptop, ratio
      REAL,   ALLOCATABLE :: deta(:), etaif(:), sfull(:), shalf(:)
      INTEGER, ALLOCATABLE :: ldum(:)
      INTEGER :: lm, lp1, k, ios, recbytes, i1, i2
!
      IF (COMMAND_ARGUMENT_COUNT() < 4) THEN
        WRITE(*,*) 'uso: gen_co2.x <deta> <PT_hPa> <ppm> <saida> [pscale|eta]'
        STOP 2
      END IF
      mode = 'pscale'
      IF (COMMAND_ARGUMENT_COUNT() >= 5) CALL GET_COMMAND_ARGUMENT(5, mode)
      CALL GET_COMMAND_ARGUMENT(1, detafile)
      CALL GET_COMMAND_ARGUMENT(2, arg); READ(arg,*) pt_hpa
      CALL GET_COMMAND_ARGUMENT(3, arg); READ(arg,*) ppm
      CALL GET_COMMAND_ARGUMENT(4, outfile)
!
!--- LM deduzido do tamanho do registro do deta: rec = 4*LM + 4*(LM+1)
      OPEN(16, FILE=detafile, FORM='UNFORMATTED', STATUS='OLD', ACTION='READ')
      INQUIRE(FILE=detafile, SIZE=recbytes)
      lm = (recbytes - 12)/8            ! arquivo = 4 + rec + 4, rec = 8*LM+4
      lp1 = lm + 1
      ALLOCATE(deta(lm), ldum(lp1), etaif(lp1), sfull(lp1), shalf(lm))
      READ(16) deta, ldum
      CLOSE(16)
      WRITE(*,'(A,I3,A,F7.2,A,F6.1,A)') ' gen_co2: LM=',lm,'  PT=',pt_hpa, &
                                        ' hPa  ppm=',ppm,'  (base 330)'
!
!--- interfaces eta (topo->fundo) e conversao p/ a convencao do SIGP
!    (fundo->topo: SFULL(1)=1 na superficie, SFULL(LP1)=0 no topo)
      etaif(1) = 0.
      DO k = 1, lm
        etaif(k+1) = etaif(k) + deta(k)
      END DO
      etaif(lp1) = 1.
      DO k = 1, lp1
        sfull(k) = etaif(lm+2-k)
      END DO
      DO k = 1, lm
        shalf(k) = 0.5*(etaif(lm+1-k) + etaif(lm+2-k))
      END DO
!
      pptop = pt_hpa/10.                ! kPa (cb), conforme SIGP/CO2PTZ
      ratio = ppm/330.
!--- modo pscale (default): concentracao via escala de pressao (RATIN=1);
!    5o argumento opcional 'eta' usa o mecanismo legado RAT->COEINT.
      IF (mode == 'eta') THEN
        CALL CO2O3(sfull, shalf, pptop, lm, lp1, lm+2, ratio, 1.0)
      ELSE
        CALL CO2O3(sfull, shalf, pptop, lm, lp1, lm+2, 1.0, ratio)
      END IF
!
!--- grava no layout do CONRAD
      OPEN(51, FILE=outfile, FORM='UNFORMATTED', STATUS='REPLACE')
      WRITE(51) (STEMP(k), k=1,lp1)
      WRITE(51) (GTEMP(k), k=1,lp1)
      WRITE(51) (CDTM51(k), k=1,lm)
      WRITE(51) (CO2M51(k), k=1,lm)
      WRITE(51) (C2DM51(k), k=1,lm)
      WRITE(51) (CDTM58(k), k=1,lm)
      WRITE(51) (CO2M58(k), k=1,lm)
      WRITE(51) (C2DM58(k), k=1,lm)
      WRITE(51) CDT51
      WRITE(51) CO251
      WRITE(51) C2D51
      WRITE(51) CDT58
      WRITE(51) CO258
      WRITE(51) C2D58
      WRITE(51) (CDT31(k), k=1,lp1)
      WRITE(51) (CO231(k), k=1,lp1)
      WRITE(51) (C2D31(k), k=1,lp1)
      WRITE(51) (CDT38(k), k=1,lp1)
      WRITE(51) (CO238(k), k=1,lp1)
      WRITE(51) (C2D38(k), k=1,lp1)
      WRITE(51) (CDT71(k), k=1,lp1)
      WRITE(51) (CO271(k), k=1,lp1)
      WRITE(51) (C2D71(k), k=1,lp1)
      WRITE(51) (CDT78(k), k=1,lp1)
      WRITE(51) (CO278(k), k=1,lp1)
      WRITE(51) (C2D78(k), k=1,lp1)
      CLOSE(51)
      WRITE(*,'(3A)') ' gen_co2: gravado ', TRIM(outfile), ' (26 registros)'
      END PROGRAM GEN_CO2
!
!--- shims p/ as chamadas wrf_* do kernel extraido
      LOGICAL FUNCTION wrf_dm_on_monitor()
      wrf_dm_on_monitor = .TRUE.
      END FUNCTION wrf_dm_on_monitor
!
      SUBROUTINE wrf_dm_bcast_bytes(buf, n)
      INTEGER :: n
      REAL    :: buf(*)
      END SUBROUTINE wrf_dm_bcast_bytes
!
      SUBROUTINE wrf_error_fatal(msg)
      CHARACTER(*) :: msg
      WRITE(*,*) 'FATAL: ', msg
      STOP 1
      END SUBROUTINE wrf_error_fatal