fix a bug in AOD calculation (#489)

* fixed a bug in AOD calculation

* add mie aod and surface dust pm10

* remove enddo in line 5211 in CLDRAD.f

* parm updates for dust pm10

parm/post_avblflds.xml

@@ -6219,16 +6219,17 @@
 
     <param>
        <post_avblfldidx>685</post_avblfldidx>
-       <shortname>DU_CR_AER_SFC_MASS_CON</shortname>
+       <shortname>DUST10_SFC_MASS_CON</shortname>
        <stats_proc>AVE</stats_proc>
        <pdstmpl>tmpl4_48</pdstmpl>
-       <pname>MASSDEN</pname>
+       <pname>PMTC</pname>
+       <fixed_sfc1_type>surface</fixed_sfc1_type>
        <aerosol_type>dust_dry</aerosol_type>
        <typ_intvl_size>smaller_than_first_limit</typ_intvl_size>
        <scale_fact_1st_size>6</scale_fact_1st_size>
        <scale_val_1st_size>10</scale_val_1st_size>
        <fixed_sfc1_type>surface</fixed_sfc1_type>
-       <scale> 9.0</scale>
+       <scale>9.0</scale>
     </param>
 
     <param>

parm/postcntrl_gefs_chem.xml

@@ -225,6 +225,12 @@
        <scale>9.0</scale>
     </param>
 
+    <param>
+       <shortname>DUST10_SFC_MASS_CON</shortname>
+       <table_info>NCEP</table_info>
+       <scale>9.0</scale>
+    </param>
+
     <param>
        <shortname>DUST25_SFC_MASS_CON</shortname>
        <table_info>NCEP</table_info>

sorc/ncep_post.fd/ALLOCATE_ALL.f

@@ -1260,7 +1260,9 @@ SUBROUTINE ALLOCATE_ALL()
         allocate(dustallcb(im,jsta_2l:jend_2u))
         allocate(ssallcb(im,jsta_2l:jend_2u))
         allocate(dustpm(im,jsta_2l:jend_2u))
+        allocate(dustpm10(im,jsta_2l:jend_2u))	
         allocate(sspm(im,jsta_2l:jend_2u))
+	allocate(maod(im,jsta_2l:jend_2u))
 !Initialization
 !$omp parallel do private(i,j)
        do j=jsta_2l,jend_2u
@@ -1295,7 +1297,9 @@ SUBROUTINE ALLOCATE_ALL()
            dustallcb(i,j)=spval
            ssallcb(i,j)=spval
            dustpm(i,j)=spval
+	   dustpm10(i,j)=spval
            sspm(i,j)=spval
+	   maod(i,j)=spval
          enddo
        enddo
       endif

sorc/ncep_post.fd/CALPW.f

@@ -273,12 +273,10 @@ SUBROUTINE CALPW(PW,IDECID)
           DO I=1,IM
             if(PINT(I,J,L+1) <spval .and. Qdum(I,J) < spval) then
              DP      = PINT(I,J,L+1) - PINT(I,J,L)
-             PW(I,J) = PW(I,J) + Qdum(I,J)*DP*GI*HTM(I,J,L)
-            IF (IDECID == 17 .or. IDECID == 20 .or. IDECID == 21) THEN
-             PW(I,J) = PW(I,J) + Qdum(I,J)*MAX(DP,0.)*GI*HTM(I,J,L)
-            ENDIF
             IF (IDECID == 19) THEN
              PW(I,J) = PW(I,J) + Qdum(I,J)
+	    ELSE
+	     PW(I,J) = PW(I,J) + Qdum(I,J)*MAX(DP,0.)*GI*HTM(I,J,L)
             ENDIF
             IF (IDECID == 14) PWS(I,J) = PWS(I,J) + QS(I,J)*DP*GI*HTM(I,J,L)
             else

sorc/ncep_post.fd/CLDRAD.f

@@ -90,7 +90,7 @@ SUBROUTINE CLDRAD
                          AVGCPRATE,                                           &
                          DUSTCB,SSCB,BCCB,OCCB,SULFCB,DUSTPM,SSPM,aod550,     &
                          du_aod550,ss_aod550,su_aod550,oc_aod550,bc_aod550,   &
-                         PWAT
+                         PWAT,DUSTPM10,MAOD
       use masks,    only: LMH, HTM
       use params_mod, only: TFRZ, D00, H99999, QCLDMIN, SMALL, D608, H1, ROG, &
                             GI, RD, QCONV, ABSCOEFI, ABSCOEF, STBOL, PQ0, A2, &
@@ -5160,6 +5160,26 @@ SUBROUTINE CLDRAD
           datapd(1:im,1:jend-jsta+1,cfld)=GRID1(1:im,jsta:jend)
          endif
       ENDIF
+      
+!! Multiply by 1.E-6 to revert these fields back
+      IF (IGET(667)>0) THEN
+         GRID1=SPVAL
+!$omp parallel do private(i,j)
+         DO J = JSTA,JEND
+            DO I = 1,IM
+               IF(BCEM(I,J,1)<SPVAL) GRID1(I,J) = BCEM(I,J,1)
+               DO K=2,NBIN_BC
+               IF(BCEM(I,J,K)<SPVAL)&
+                GRID1(I,J) = GRID1(I,J) + BCEM(I,J,K)
+               END DO
+            END DO
+         END DO
+         if(grib=='grib2') then
+          cfld=cfld+1
+          fld_info(cfld)%ifld=IAVBLFLD(IGET(667))
+          datapd(1:im,1:jend-jsta+1,cfld)=GRID1(1:im,jsta:jend)
+         endif
+      ENDIF
 
       IF (IGET(660)>0) THEN
          GRID1=SPVAL
@@ -5179,6 +5199,21 @@ SUBROUTINE CLDRAD
           datapd(1:im,1:jend-jsta+1,cfld)=GRID1(1:im,jsta:jend)
          endif
       ENDIF
+      
+      IF (IGET(699)>0) THEN
+         GRID1=SPVAL
+!$omp parallel do private(i,j)
+         DO J = JSTA,JEND
+            DO I = 1,IM
+                GRID1(I,J) = MAOD(I,J)
+            END DO
+         END DO
+         if(grib=='grib2') then
+          cfld=cfld+1
+          fld_info(cfld)%ifld=IAVBLFLD(IGET(699))
+          datapd(1:im,1:jend-jsta+1,cfld)=GRID1(1:im,jsta:jend)
+         endif
+      ENDIF
 !! ADD DUST DRY DEPOSITION FLUXES (kg/m2/sec)
 !
 !      IF (IGET(661)>0) THEN
@@ -5217,6 +5252,20 @@ SUBROUTINE CLDRAD
          endif
       ENDIF
 
+      IF (IGET(685)>0 ) THEN
+!$omp parallel do private(i,j)
+         DO J = JSTA,JEND
+            DO I = 1,IM
+               GRID1(I,J) = DUSTPM10(I,J)   !ug/m3 
+            END DO
+         END DO
+         if(grib=='grib2') then
+           cfld=cfld+1
+           fld_info(cfld)%ifld=IAVBLFLD(IGET(685))
+           datapd(1:im,1:jend-jsta+1,cfld)=GRID1(1:im,jsta:jend)
+         endif
+      ENDIF
+      
 !! ADD DUST WET DEPOSITION FLUXES (kg/m2/sec)
 !      IF (IGET(662)>0) THEN
 !         DO J = JSTA,JEND
@@ -5430,7 +5479,10 @@ SUBROUTINE CLDRAD
       IF (IGET(674)>0) call wrt_aero_diag(674,nbin_oc,ocwt)
       IF (IGET(682)>0) call wrt_aero_diag(682,nbin_oc,ocsv)
 !      print *,'aft wrt disg ocwt'
-
+!! wrt MIE AOD at 550nm
+      IF (IGET(699).GT.0) call wrt_aero_diag(699,1,maod)
+      print *,'aft wrt disg maod'
+      
 !! wrt SU diag field
 !      IF (IGET(675)>0) call wrt_aero_diag(675,nbin_su,suem)
 !      IF (IGET(676)>0) call wrt_aero_diag(676,nbin_su,susd)

sorc/ncep_post.fd/DEALLOCATE.f

@@ -532,7 +532,9 @@ SUBROUTINE DE_ALLOCATE
         deallocate(dustallcb)
         deallocate(ssallcb)
         deallocate(dustpm)
+	deallocate(dustpm10)
         deallocate(sspm)
+	deallocate(maod)
       endif
 !
 ! HWRF RRTMG output