Changeset 11d86e7 in flexpart.git
- Timestamp:
- May 28, 2020, 7:16:49 AM (4 years ago)
- Branches:
- dev
- Children:
- 969c1fe
- Parents:
- ffbe224
- Location:
- src
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
src/netcdf_output_mod.f90
rffbe224 r11d86e7 33 33 use outg_mod, only: outheight,oroout,densityoutgrid,factor3d,volume,& 34 34 wetgrid,wetgridsigma,drygrid,drygridsigma,grid,gridsigma,& 35 area,arean,volumen, orooutn, areaeast 35 area,arean,volumen, orooutn, areaeast, areanorth 36 36 use par_mod, only: dep_prec, sp, dp, maxspec, maxreceptor, nclassunc,& 37 37 unitoutrecept,unitoutreceptppt, nxmax,unittmp … … 1566 1566 call nf90_err(nf90_put_var(ncid,fluxid,1.e12*grid(0:numxgrid-1,0:numygrid-1,1:numzgrid)& 1567 1567 /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 /) )) 1569 1649 1570 1650 end do -
src/timemanager.f90
rffbe224 r11d86e7 430 430 if ((iout.eq.4).or.(iout.eq.5)) call plumetraj(itime) 431 431 if ((iflux.eq.1).and.(lnetcdfout.eq.0)) call fluxoutput(itime) 432 #ifdef USE_NCF 432 433 if ((iflux.eq.1).and.(lnetcdfout.eq.1)) call fluxoutput_netcdf(itime) 434 #endif 433 435 write(*,45) itime,numpart,gridtotalunc,wetgridtotalunc,drygridtotalunc 434 436
Note: See TracChangeset
for help on using the changeset viewer.