Changeset 11d86e7 in flexpart.git


Ignore:
Timestamp:
May 28, 2020, 7:16:49 AM (4 years ago)
Author:
Sabine <sabine.eckhardt@…>
Branches:
dev
Children:
969c1fe
Parents:
ffbe224
Message:

bugfix in the fluxoutput

Location:
src
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • src/netcdf_output_mod.f90

    rffbe224 r11d86e7  
    3333  use outg_mod,  only: outheight,oroout,densityoutgrid,factor3d,volume,&
    3434                       wetgrid,wetgridsigma,drygrid,drygridsigma,grid,gridsigma,&
    35                        area,arean,volumen, orooutn, areaeast
     35                       area,arean,volumen, orooutn, areaeast, areanorth
    3636  use par_mod,   only: dep_prec, sp, dp, maxspec, maxreceptor, nclassunc,&
    3737                       unitoutrecept,unitoutreceptppt, nxmax,unittmp
     
    15661566        call nf90_err(nf90_put_var(ncid,fluxid,1.e12*grid(0:numxgrid-1,0:numygrid-1,1:numzgrid)&
    15671567            /areaeast(0:numxgrid-1,0:numygrid-1,1:numzgrid)/loutstep,&
    1568            (/ 1,1,1,tpointer,kp,nage /), (/ numxgrid,numygrid,numzgrid,1,1,1 /) ))
     1568           (/ 1,1,1,1,kp,nage /), (/ numxgrid,numygrid,numzgrid,1,1,1 /) ))
     1569
     1570 ! West Flux
     1571        call nf90_err(nf90_def_var(ncid,'flux_west_'//anspec, nf90_float, dIDs, &
     1572             fluxID))
     1573
     1574        do jy=0,numygrid-1
     1575          do ix=0,numxgrid-1
     1576            do kz=1, numzgrid
     1577               grid(ix,jy,kz)=flux(2,ix,jy,kz,ks,kp,nage)
     1578            end do
     1579          end do
     1580        end do
     1581
     1582        call nf90_err(nf90_put_var(ncid,fluxid,1.e12*grid(0:numxgrid-1,0:numygrid-1,1:numzgrid)&
     1583            /areaeast(0:numxgrid-1,0:numygrid-1,1:numzgrid)/loutstep,&
     1584           (/ 1,1,1,1,kp,nage /), (/ numxgrid,numygrid,numzgrid,1,1,1 /) ))
     1585
     1586 ! North Flux
     1587        call nf90_err(nf90_def_var(ncid,'flux_north_'//anspec, nf90_float, dIDs, &
     1588             fluxID))
     1589
     1590        do jy=0,numygrid-1
     1591          do ix=0,numxgrid-1
     1592            do kz=1, numzgrid
     1593               grid(ix,jy,kz)=flux(4,ix,jy,kz,ks,kp,nage)
     1594            end do
     1595          end do
     1596        end do
     1597
     1598        call nf90_err(nf90_put_var(ncid,fluxid,1.e12*grid(0:numxgrid-1,0:numygrid-1,1:numzgrid)&
     1599            /areanorth(0:numxgrid-1,0:numygrid-1,1:numzgrid)/loutstep,&
     1600           (/ 1,1,1,1,kp,nage /), (/ numxgrid,numygrid,numzgrid,1,1,1 /) ))
     1601
     1602 ! South Flux
     1603        call nf90_err(nf90_def_var(ncid,'flux_south_'//anspec, nf90_float, dIDs, &
     1604             fluxID))
     1605
     1606        do jy=0,numygrid-1
     1607          do ix=0,numxgrid-1
     1608            do kz=1, numzgrid
     1609               grid(ix,jy,kz)=flux(3,ix,jy,kz,ks,kp,nage)
     1610            end do
     1611          end do
     1612        end do
     1613
     1614        call nf90_err(nf90_put_var(ncid,fluxid,1.e12*grid(0:numxgrid-1,0:numygrid-1,1:numzgrid)&
     1615            /areanorth(0:numxgrid-1,0:numygrid-1,1:numzgrid)/loutstep,&
     1616           (/ 1,1,1,1,kp,nage /), (/ numxgrid,numygrid,numzgrid,1,1,1 /) ))
     1617
     1618 ! Up Flux
     1619        call nf90_err(nf90_def_var(ncid,'flux_up_'//anspec, nf90_float, dIDs, &
     1620             fluxID))
     1621
     1622        do jy=0,numygrid-1
     1623          do ix=0,numxgrid-1
     1624            do kz=1, numzgrid
     1625               grid(ix,jy,kz)=flux(5,ix,jy,kz,ks,kp,nage)/area(ix,jy)
     1626            end do
     1627          end do
     1628        end do
     1629
     1630        call nf90_err(nf90_put_var(ncid,fluxid,1.e12*grid(0:numxgrid-1,0:numygrid-1,1:numzgrid)&
     1631            /loutstep,&
     1632           (/ 1,1,1,1,kp,nage /), (/ numxgrid,numygrid,numzgrid,1,1,1 /) ))
     1633
     1634 ! Down Flux
     1635        call nf90_err(nf90_def_var(ncid,'flux_down_'//anspec, nf90_float, dIDs, &
     1636             fluxID))
     1637
     1638        do jy=0,numygrid-1
     1639          do ix=0,numxgrid-1
     1640            do kz=1, numzgrid
     1641               grid(ix,jy,kz)=flux(6,ix,jy,kz,ks,kp,nage)/area(ix,jy)
     1642            end do
     1643          end do
     1644        end do
     1645
     1646        call nf90_err(nf90_put_var(ncid,fluxid,1.e12*grid(0:numxgrid-1,0:numygrid-1,1:numzgrid)&
     1647            /loutstep,&
     1648           (/ 1,1,1,1,kp,nage /), (/ numxgrid,numygrid,numzgrid,1,1,1 /) ))
    15691649
    15701650      end do
  • src/timemanager.f90

    rffbe224 r11d86e7  
    430430        if ((iout.eq.4).or.(iout.eq.5)) call plumetraj(itime)
    431431        if ((iflux.eq.1).and.(lnetcdfout.eq.0)) call fluxoutput(itime)
     432#ifdef USE_NCF
    432433        if ((iflux.eq.1).and.(lnetcdfout.eq.1)) call fluxoutput_netcdf(itime)
     434#endif
    433435        write(*,45) itime,numpart,gridtotalunc,wetgridtotalunc,drygridtotalunc
    434436 
Note: See TracChangeset for help on using the changeset viewer.
hosted by ZAMG