      SUBROUTINE NM2PLA (PN0,PHI0,PSI,RESN,RESM,BNUM,BMASS, &
                         DRYDN,PHISS,DPHIS,PHIS0,DPHI0, &
                         ILGA,LEVA,ISEC,KPT)
!-----------------------------------------------------------------------
!     PURPOSE:
!     --------
!     TRANSFORMATION OF NUMBER AND MASS TO PLA-PARAMETERS. SUBROUTINE
!     CORNM SHOULD BE CALLED BEFORE THIS SUBROUTINE.
!
!     HISTORY:
!     --------
!     * JAN 07/2008 - K.VONSALZEN   IMPROVE NUMERICAL STABILITY.
!     * APR 13/2006 - K.VONSALZEN   COMPLETELY REVISED FOR USE WITH
!     *                             INTERNALLY AND EXTERNALLY MIXED
!     *                             TYPES OF AEROSOLS.
!     * DEC 11/2005 - K.VONSALZEN   NEW.
!
!     ************************ INDEX OF VARIABLES **********************
!
!     PRIMARY INPUT/OUTPUT VARIABLES (ARGUMENT LIST):
!
! I   BMASS    AEROSOL MASS CONCENTRATION (KG/KG)
! I   BNUM     AEROSOL NUMBER CONCENTRATION (1/KG)
! I   DPHI0    SECTION WIDTH
! I   DPHIS    SECTION WIDTH
! I   DRYDN    DENSITY OF DRY AEROSOL PARTICLE (KG/M3)
! I   ILGA     NUMBER OF GRID POINTS IN HORIZONTAL DIRECTION
! I   ISEC     NUMBER OF SEPARATE AEROSOL TRACERS
! I   KPT      FLAG TO DISTINGUISH BETWEEN EXTERNALLY
!              RESP. INTERNALLY MIXED TYPES OF AEROSOL
! I   LEVA     NUMBER OF GRID POINTA IN VERTICAL DIRECTION
! O   PHI0     3RD PLA SIZE DISTRIBUTION PARAMETER (MODE SIZE)
! I   PHIS0    LOWER DRY PARTICLE SIZE (=LOG(R/R0))
! I   PHISS    LOWER DRY PARTICLE SIZE (=LOG(R/R0))
! O   PN0      1ST PLA SIZE DISTRIBUTION PARAMETER (AMPLITUDE)
! O   PSI      2ND PLA SIZE DISTRIBUTION PARAMETER (WIDTH)
! I   PTRANS(:)%FLG(:)
!              FLAG TO INDICATE VALID DATA IN PLA DATA TABLE
! I   PTRANS(:)%IDRB
!              NUMBER OF POINTS IN PLA DATA TABLE
! I   PTRANS(:)%DPHIB(:)%VL
!              TABULATED DERIVATIVE OF DELTA-PHI AT LOWER SECTION
!              BOUNDARY
! I   PTRANS(:)%DPHIB(:)%VR
!              TABULATED DERIVATIVE OF DELTA-PHI AT UPPER SECTION
!              BOUNDARY
! I   PTRANS(:)%PHIB(:)%VL
!              TABULATED RESULTS FOR DELTA-PHI AT LOWER SECTION
!              BOUNDARY
! I   PTRANS(:)%PHIB(:)%VR
!              TABULATED RESULTS FOR DELTA-PHI AT UPPER SECTION
!              BOUNDARY
! I   PTRANS(:)%PSIB(:)%VL
!              TABULATED RESULTS FOR PSI AT LOWER SECTION BOUNDARY
! O   RESM     MASS RESIDUUM FROM NUMERICAL TRUNCATION (KG/KG)
! O   RESN     NUMBER RESIDUUM FROM NUMERICAL TRUNCATION (KG/KG)
!
!     SECONDARY INPUT/OUTPUT VARIABLES (VIA MODULES AND COMMON BLOCKS):
!
! I   ILARGE   INTEGER EQUIVALENT TO YLARGE
! I   R0       PARTICLE RADIUS REFERENCE RADIUS (1.E-06 M)
! I   YCNST    4*PI/3
! I   YLARGE   LARGE NUMBER
! I   YNA      DEFAULT FOR UNDEFINED VALUE
! I   YTINY    SMALLEST VALID NUMBER
!
!-----------------------------------------------------------------------
!
      USE SDPARM
      USE SDCODE
      USE FPDEF
!
      IMPLICIT REAL (A-H,O-Z), INTEGER (I-N)
!
      LOGICAL, INTENT(IN) :: KPT
      INTEGER, INTENT(IN) :: ILGA,LEVA,ISEC
      REAL, INTENT(IN), DIMENSION(ISEC) :: PHISS,DPHIS
      REAL, INTENT(IN), DIMENSION(ILGA,LEVA,ISEC) :: BNUM,BMASS,DRYDN, &
                                                     PHIS0,DPHI0
      REAL, INTENT(OUT), DIMENSION(ILGA,LEVA) :: RESN,RESM
      REAL, INTENT(OUT), DIMENSION(ILGA,LEVA,ISEC) :: PN0,PHI0,PSI
      REAL, ALLOCATABLE, DIMENSION(:,:,:) :: PN0N,PN0M,SDFN,SDFM,RAT0S, &
                                             ARG,PHIHAT,ATMP,AMASS
      REAL, ALLOCATABLE, DIMENSION(:,:) :: TOTN,TOTM
      TYPE(PTRANT), ALLOCATABLE, DIMENSION(:) :: PTRANS
!
      REAL, PARAMETER :: ONETHIRD=1./3.
      REAL, PARAMETER :: YSCAL=1.E+09
      REAL, PARAMETER :: RWEIGH=1.         ! WEIGHTING FACTOR
      INTEGER, PARAMETER :: ITMAX=5        ! NUMBER OF ITERATIONS FOR
                                           ! INCREASED ACCURACY
      INTEGER, PARAMETER :: ICONSM=2       ! DETERMINES IF NUMERICAL
                                           ! TRUNCATION WILL PRIMARILY
                                           ! AFFECT MASS OR NUMBER RESULTS
      LOGICAL, PARAMETER :: KITER=.FALSE.  ! ITERATIONS FOR INCREASED
                                           ! ACCURACY
      REAL(R8) :: RMOM3,TRM1,TRM2,SQRTPSI,DPHI,ARG1,ARG2,ARG3,ARG4, &
                  AERF1,AERF2,AERF3,AERF4,ATERM1,ADERF1,ADERF2, &
                  ADERF3,ADERF4,ADFUN,AFUN
!
!-----------------------------------------------------------------------
!     * ALLOCATE WORK ARRAYS.
!
      ALLOCATE(PN0N  (ILGA,LEVA,ISEC))
      ALLOCATE(PN0M  (ILGA,LEVA,ISEC))
      ALLOCATE(SDFN  (ILGA,LEVA,ISEC))
      ALLOCATE(SDFM  (ILGA,LEVA,ISEC))
      ALLOCATE(RAT0S (ILGA,LEVA,ISEC))
      ALLOCATE(ARG   (ILGA,LEVA,ISEC))
      ALLOCATE(PHIHAT(ILGA,LEVA,ISEC))
      ALLOCATE(ATMP  (ILGA,LEVA,ISEC))
      ALLOCATE(AMASS (ILGA,LEVA,ISEC))
      ALLOCATE(PTRANS(ISATR))
      DO IS=1,ISATR
        ALLOCATE(PTRANS(IS)%FLG  (IDRBM))
        ALLOCATE(PTRANS(IS)%PSIB (IDRBM))
        ALLOCATE(PTRANS(IS)%PHIB (IDRBM))
        ALLOCATE(PTRANS(IS)%DPHIB(IDRBM))
      ENDDO
!
!-----------------------------------------------------------------------
!     * INITIALIZATIONS.
!
      PN0 =0.
      PN0M=0.
      PN0N=0.
      PSI =YNA
      PHI0=YNA
      R0P=R0**3
      ACNST=YCNST*R0P*YSCAL
!
!     * COPY PARAMETERS SPECIFIC TO INTERNALLY RESP. EXTERNALLY
!     * MIXED TYPES OF AEROSOL FROM PERMANENT TO TEMPORARY MEMORY.
!
      IF ( KPT .EQV. KPEXT ) THEN
        DO IS=1,ISATR
          PTRANS(IS)%IDRB       =PEXT(IS)%IDRB
          PTRANS(IS)%FLG  (:)   =PEXT(IS)%FLG  (:)
          PTRANS(IS)%PSIB (:)%VL=PEXT(IS)%PSIB (:)%VL
          PTRANS(IS)%PSIB (:)%VR=PEXT(IS)%PSIB (:)%VR
          PTRANS(IS)%PHIB (:)%VL=PEXT(IS)%PHIB (:)%VL
          PTRANS(IS)%PHIB (:)%VR=PEXT(IS)%PHIB (:)%VR
          PTRANS(IS)%DPHIB(:)%VL=PEXT(IS)%DPHIB(:)%VL
          PTRANS(IS)%DPHIB(:)%VR=PEXT(IS)%DPHIB(:)%VR
        ENDDO
      ELSE
        DO IS=1,ISATR
          PTRANS(IS)%IDRB       =PINT(IS)%IDRB
          PTRANS(IS)%FLG  (:)   =PINT(IS)%FLG  (:)
          PTRANS(IS)%PSIB (:)%VL=PINT(IS)%PSIB (:)%VL
          PTRANS(IS)%PSIB (:)%VR=PINT(IS)%PSIB (:)%VR
          PTRANS(IS)%PHIB (:)%VL=PINT(IS)%PHIB (:)%VL
          PTRANS(IS)%PHIB (:)%VR=PINT(IS)%PHIB (:)%VR
          PTRANS(IS)%DPHIB(:)%VL=PINT(IS)%DPHIB(:)%VL
          PTRANS(IS)%DPHIB(:)%VR=PINT(IS)%DPHIB(:)%VR
        ENDDO
      ENDIF
!
!-----------------------------------------------------------------------
!     * VALUES OF PSI AND PHI0.
!
      ARG=1.
      ATMP=ACNST*DRYDN*BNUM
      AMASS=BMASS*YSCAL
      WHERE ( BNUM > YTINY .AND. BMASS > YTINY &
                                 .AND. ATMP > MAX(AMASS/YLARGE,YTINY) )
        ARG=AMASS/ATMP
      ENDWHERE
      PHIHAT=ONETHIRD*LOG(ARG)
      DO IS=1,ISEC
        IDRB=PTRANS(IS)%IDRB
        DRATB=1./REAL(IDRB)
        DO L=1,LEVA
        DO IL=1,ILGA
          IF ( BNUM(IL,L,IS) > YTINY .AND. BMASS(IL,L,IS) > YTINY &
                                                                ) THEN
!
!           * SIZE RATIO FROM MASS AND NUMBER IN EACH SECTION.
!
            RAT0=(PHIHAT(IL,L,IS)-PHISS(IS))/DPHIS(IS)
            RAT0S(IL,L,IS)=MAX(0.,MIN(1.,RAT0))
            IR=MIN(NINT(RAT0*IDRB+0.5),IDRB)
!
!           * CHECK IF SIZE RATIO IS VALID. ABORT IF THIS IS NOT THE CASE
!           * SUBROUTINE CORNM WILL PREVENT THIS FROM HAPPENING IF
!           * CALLED BEFORE THIS SUBROUTINE.
!
            IF ( INT(RAT0) /= 0 .OR. RAT0 < 0. ) CALL XIT ('NM2PLA',-1)
            IF ( PTRANS(IS)%FLG(IR) ) THEN
!
!             * USE TABULATED RESULTS TO CALCULATE PHI-PARAMETER FROM
!             * LINEAR INTERPOLATION AND TO SET PSI-PARAMETER.
!
              PSI(IL,L,IS)=PTRANS(IS)%PSIB(IR)%VL
              DRATX=RAT0-REAL(IR-1)*DRATB
              WEIGF=1.-DRATX/DRATB
              PHILX=PTRANS(IS)%PHIB(IR)%VL &
                                        +DRATX*PTRANS(IS)%DPHIB(IR)%VL
              PHIRX=PTRANS(IS)%PHIB(IR)%VR &
                                -(DRATB-DRATX)*PTRANS(IS)%DPHIB(IR)%VR
              PHI0(IL,L,IS)=PHISS(IS)+WEIGF*PHILX+(1.-WEIGF)*PHIRX
            ELSE
!
!             * SIZE DISTRIBUTION IS A DELTA-FUNCTION
!
              PSI(IL,L,IS)=YLARGE
              PHI0(IL,L,IS)=PHIHAT(IL,L,IS)
            ENDIF
          ENDIF
        ENDDO
        ENDDO
      ENDDO
!
!-----------------------------------------------------------------------
!     * DEALLOCATE WORK ARRAYS.
!
      DEALLOCATE(ARG)
      DEALLOCATE(PHIHAT)
      DEALLOCATE(PTRANS)
!
!-----------------------------------------------------------------------
!     * ITERATIONS TO INCREASE ACCURACY.
!
      IF ( KITER ) THEN
        RMOM=3.
        TRM1=.5*RMOM
        TRM2=0.25*RMOM**2
        DO IS=1,ISEC
        DO L=1,LEVA
        DO IL=1,ILGA
          IF ( PSI(IL,L,IS) < YLARGES &
                             .AND. ABS(PSI(IL,L,IS)-YNA) > YTINY ) THEN
            SQRTPSI=SQRT(ABS(PSI(IL,L,IS)))
            DPHI=PHI0(IL,L,IS)-PHISS(IS)
            IOK=1
            DO IT=1,ITMAX
              ARG3= SQRTPSI*(DPHIS(IS)-DPHI)
              ARG4=-SQRTPSI*DPHI
              IF ( PSI(IL,L,IS) > YTINY ) THEN
                ARG1=ARG3-TRM1/SQRTPSI
                ARG2=ARG4-TRM1/SQRTPSI
                AERF1=                                                 &!PROC-DEPEND
                      ERF(ARG1)
                AERF2=                                                 &!PROC-DEPEND
                      ERF(ARG2)
                AERF3=                                                 &!PROC-DEPEND
                      ERF(ARG3)
                AERF4=                                                 &!PROC-DEPEND
                      ERF(ARG4)
                SIGP=1.
              ELSE
                ARG1=ARG3+TRM1/SQRTPSI
                ARG2=ARG4+TRM1/SQRTPSI
                AERF1=ERFI(ARG1)
                AERF2=ERFI(ARG2)
                AERF3=ERFI(ARG3)
                AERF4=ERFI(ARG4)
                SIGP=-1.
              ENDIF
!
!             * CRITICAL FUNCTION.
!
              ARG1=3.*(RAT0S(IL,L,IS)*DPHIS(IS)-DPHI) &
                                                  -9./(4.*PSI(IL,L,IS))
              ATERM1=                                                  &!PROC-DEPEND
                     EXP(ARG1)
              IF ( ABS(AERF3-AERF4) > YTINY ) THEN
                ATERM2=(AERF1-AERF2)/(AERF3-AERF4)
                AFUN=ATERM1-ATERM2
!
!               * DERIVATIVE.
!
                ARG1=-SIGP*ARG1**2
                ARG2=-SIGP*ARG2**2
                ARG3=-SIGP*ARG3**2
                ARG4=-SIGP*ARG4**2
                ADERF1=-2.*                                            &!PROC-DEPEND
                           EXP(ARG1)*SQRTPSI/SQRTPI
                ADERF2=-2.*                                            &!PROC-DEPEND
                           EXP(ARG2)*SQRTPSI/SQRTPI
                ADERF3=-2.*                                            &!PROC-DEPEND
                           EXP(ARG3)*SQRTPSI/SQRTPI
                ADERF4=-2.*                                            &!PROC-DEPEND
                           EXP(ARG4)*SQRTPSI/SQRTPI
                ADFUN=-3.*ATERM1-((ADERF1-ADERF2)*(AERF3-AERF4) &
                                 -(ADERF3-ADERF4)*(AERF1-AERF2)) &
                                /(AERF3-AERF4)**2
!
!               * UPDATE DPHI
!
                DPHI=DPHI-AFUN/ADFUN
              ELSE
                IOK=IDEF
              ENDIF
            ENDDO
            IF ( IOK /= IDEF ) PHI0(IL,L,IS)=PHISS(IS)+DPHI
          ENDIF
        ENDDO
        ENDDO
        ENDDO
      ENDIF
!
!-----------------------------------------------------------------------
!     * DEALLOCATE WORK ARRAY.
!
      DEALLOCATE(RAT0S)
!
!-----------------------------------------------------------------------
!     * CALCULATE REMAINING THIRD PARAMETER OF GAUSSIAN DISTRIBUTION.
!     * THE CALCULATION OF THE PHI-PARAMETER FROM TABULATED RESULTS
!     * WILL LEAD TO INACCURACIES FOR THIS PARAMETER.
!
      SDFN=SDINT0(PHI0,PSI,PHIS0,DPHI0,ILGA,LEVA,ISEC)
      RMOM=3.
      SDFM=DRYDN*YCNST*SDINT(PHI0,PSI,RMOM,PHIS0,DPHI0,ILGA,LEVA,ISEC)
      WHERE (          SDFN > MAX(BNUM /YLARGE,YTINY) &
                 .AND. SDFM > MAX(BMASS/YLARGE,YTINY) )
        PN0N=BNUM /SDFN
        PN0M=BMASS/SDFM
      ENDWHERE
      IF ( ICONSM==1 ) THEN
!
!       * TRUNCATION ERRORS ARE ASSUMED TO ONLY AFFECT THE MASS
!       * SIZE DISTRIBUTION.
!
        PN0=PN0N
        RESN=0.
        RESM=SUM(PN0*SDFM-BMASS,DIM=3)
      ELSE IF ( ICONSM==2 ) THEN
!
!       * TRUNCATION ERRORS ARE ASSUMED TO ONLY AFFECT THE NUMBER
!       * SIZE DISTRIBUTION.
!
        PN0=PN0M
        RESN=SUM(PN0*SDFN-BNUM ,DIM=3)
        RESM=0.
      ELSE IF ( ICONSM==3 ) THEN
!
!       * REDUCE THE IMPACT OF INACCURACIES ON TOTAL MASS AND NUMBER
!       * BY USING A WEIGHTED AVERAGE OF THE PARAMETERS THAT WOULD
!       * RESULT FOR MASS AND NUMBER ONLY.
!
        ALLOCATE(TOTN(ILGA,LEVA))
        ALLOCATE(TOTM(ILGA,LEVA))
        TOTM=SUM(BMASS,DIM=3)
        TOTN=SUM(BNUM ,DIM=3)
        DO IS=1,ISEC
        DO L=1,LEVA
        DO IL=1,ILGA
          IF ( BNUM(IL,L,IS) > YTINY .AND. BMASS(IL,L,IS) > YTINY &
                                                                ) THEN
            FRCM=RWEIGH*BMASS(IL,L,IS)/TOTM(IL,L)
            FRCN=BNUM(IL,L,IS)/TOTN(IL,L)
            WEIGF=FRCN/(FRCN+FRCM)
            PN0(IL,L,IS)=WEIGF*PN0N(IL,L,IS)+(1.-WEIGF)*PN0M(IL,L,IS)
          ENDIF
        ENDDO
        ENDDO
        ENDDO
        RESN=SUM(PN0*SDFN-BNUM ,DIM=3)
        RESM=SUM(PN0*SDFM-BMASS,DIM=3)
        DEALLOCATE(TOTN)
        DEALLOCATE(TOTM)
      ENDIF
!
!-----------------------------------------------------------------------
!     * DEALLOCATE WORK ARRAYS.
!
      DEALLOCATE(PN0N )
      DEALLOCATE(PN0M )
      DEALLOCATE(SDFN )
      DEALLOCATE(SDFM )
      DEALLOCATE(ATMP )
      DEALLOCATE(AMASS)
!
      END SUBROUTINE NM2PLA
