Changeset ffbe224 in flexpart.git
- Timestamp:
- Jan 20, 2020, 10:01:44 AM (4 years ago)
- Branches:
- dev
- Children:
- 11d86e7
- Parents:
- 95a8cb6
- Location:
- src
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
src/init_domainfill.f90
r3481cc1 rffbe224 206 206 ztra1(numpart+jj)=height(nz)-0.5 207 207 208 209 208 ! Interpolate PV to the particle position 210 209 !**************************************** … … 213 212 ixp=ixm+1 214 213 jyp=jym+1 214 if (jyp.gt.180) then 215 write (*,*) 'init_domainfill, over: ',jyp,jym,ytra1(numpart+jj),jy,ran1(idummy),ny 216 jyp=jym 217 endif 215 218 ddx=xtra1(numpart+jj)-real(ixm) 216 219 ddy=ytra1(numpart+jj)-real(jym) -
src/netcdf_output_mod.f90
r3481cc1 rffbe224 33 33 use outg_mod, only: outheight,oroout,densityoutgrid,factor3d,volume,& 34 34 wetgrid,wetgridsigma,drygrid,drygridsigma,grid,gridsigma,& 35 area,arean,volumen, orooutn 35 area,arean,volumen, orooutn, areaeast 36 36 use par_mod, only: dep_prec, sp, dp, maxspec, maxreceptor, nclassunc,& 37 37 unitoutrecept,unitoutreceptppt, nxmax,unittmp 38 use com_mod, only: path,length,ldirect, ibdate,ibtime,iedate,ietime, &38 use com_mod, only: path,length,ldirect,bdate,ibdate,ibtime,iedate,ietime, & 39 39 loutstep,loutaver,loutsample,outlon0,outlat0,& 40 40 numxgrid,numygrid,dxout,dyout,numzgrid, height, & … … 63 63 64 64 public :: writeheader_netcdf,concoutput_surf_nest_netcdf,concoutput_netcdf,& 65 &concoutput_nest_netcdf,concoutput_surf_netcdf 65 &concoutput_nest_netcdf,concoutput_surf_netcdf,fluxoutput_netcdf 66 66 67 67 ! include 'netcdf.inc' … … 1424 1424 end subroutine concoutput_surf_nest_netcdf 1425 1425 1426 subroutine fluxoutput_netcdf(itime) 1427 1428 ! i 1429 !***************************************************************************** 1430 ! * 1431 ! Output of the gridded fluxes. * 1432 ! Eastward, westward, northward, southward, upward and downward gross * 1433 ! fluxes are written to output file in either sparse matrix or grid dump * 1434 ! format, whichever is more efficient. * 1435 ! * 1436 ! Author: A. Stohl * 1437 ! * 1438 ! 04 April 2000 * 1439 ! netcdfoutput S. Eckhardt, 2020 * 1440 ! * 1441 !***************************************************************************** 1442 1443 use flux_mod 1444 1445 implicit none 1446 1447 real(kind=dp) :: jul 1448 integer :: itime,ix,jy,kz,ks,nage,jjjjmmdd,ihmmss,kp,i 1449 1450 character(len=255) :: ncfname 1451 character :: adate*8,atime*6,timeunit*32,anspec*3 1452 1453 integer :: ncid 1454 integer :: timeDimID, latDimID, lonDimID, levDimID 1455 integer :: nspecDimID, npointDimID, nageclassDimID, ncharDimID, pointspecDimID 1456 integer :: tID, lonID, latID, levID, lageID, fluxID 1457 integer, dimension(6) :: dIDs 1458 integer :: cache_size 1459 real, allocatable, dimension(:) :: coord 1460 1461 1462 ! Determine current calendar date, needed for the file name 1463 !********************************************************** 1464 1465 jul=bdate+real(itime,kind=dp)/86400._dp 1466 call caldate(jul,jjjjmmdd,ihmmss) 1467 write(adate,'(i8.8)') jjjjmmdd 1468 write(atime,'(i6.6)') ihmmss 1469 1470 ncfname=path(2)(1:length(2))//'grid_flux_'//adate// & 1471 atime//'.nc' 1472 1473 ! setting cache size in bytes. It is set to 4 times the largest data block that is written 1474 ! size_type x nx x ny x nz 1475 ! create file 1476 1477 cache_size = 16 * numxgrid * numygrid * numzgrid 1478 call nf90_err(nf90_create(trim(ncfname), cmode = nf90_hdf5, ncid = ncid, & 1479 cache_size = cache_size)) 1480 1481 ! create dimensions: 1482 !************************* 1483 ! time 1484 call nf90_err(nf90_def_dim(ncid, 'time', nf90_unlimited, timeDimID)) 1485 timeunit = 'seconds since '//adate(1:4)//'-'//adate(5:6)// & 1486 '-'//adate(7:8)//' '//atime(1:2)//':'//atime(3:4) 1487 1488 call nf90_err(nf90_def_dim(ncid, 'longitude', numxgrid, lonDimID)) 1489 call nf90_err(nf90_def_dim(ncid, 'latitude', numygrid, latDimID)) 1490 call nf90_err(nf90_def_dim(ncid, 'height', numzgrid, levDimID)) 1491 call nf90_err(nf90_def_dim(ncid, 'numspec', nspec, nspecDimID)) 1492 call nf90_err(nf90_def_dim(ncid, 'pointspec', maxpointspec_act, pointspecDimID)) 1493 call nf90_err(nf90_def_dim(ncid, 'nageclass', nageclass, nageclassDimID)) 1494 call nf90_err(nf90_def_dim(ncid, 'nchar', 45, ncharDimID)) 1495 call nf90_err(nf90_def_dim(ncid, 'numpoint', numpoint, npointDimID)) 1496 1497 1498 ! create variables 1499 !************************* 1500 1501 ! time 1502 call nf90_err(nf90_def_var(ncid, 'time', nf90_int, (/ timeDimID /), tID)) 1503 call nf90_err(nf90_put_att(ncid, tID, 'units', timeunit)) 1504 call nf90_err(nf90_put_att(ncid, tID, 'calendar', 'proleptic_gregorian')) 1505 timeID = tID 1506 1507 ! lon 1508 call nf90_err(nf90_def_var(ncid, 'longitude', nf90_float, (/ lonDimID /), lonID)) 1509 call nf90_err(nf90_put_att(ncid, lonID, 'long_name', 'longitude in degree east')) 1510 call nf90_err(nf90_put_att(ncid, lonID, 'axis', 'Lon')) 1511 call nf90_err(nf90_put_att(ncid, lonID, 'units', 'degrees_east')) 1512 call nf90_err(nf90_put_att(ncid, lonID, 'standard_name', 'grid_longitude')) 1513 call nf90_err(nf90_put_att(ncid, lonID, 'description', 'grid cell centers')) 1514 1515 ! lat 1516 call nf90_err(nf90_def_var(ncid, 'latitude', nf90_float, (/ latDimID /), latID)) 1517 call nf90_err(nf90_put_att(ncid, latID, 'long_name', 'latitude in degree north')) 1518 call nf90_err(nf90_put_att(ncid, latID, 'axis', 'Lat')) 1519 call nf90_err(nf90_put_att(ncid, latID, 'units', 'degrees_north')) 1520 call nf90_err(nf90_put_att(ncid, latID, 'standard_name', 'grid_latitude')) 1521 call nf90_err(nf90_put_att(ncid, latID, 'description', 'grid cell centers')) 1522 1523 ! height 1524 call nf90_err(nf90_def_var(ncid, 'height', nf90_float, (/ levDimID /), levID)) 1525 call nf90_err(nf90_put_att(ncid, levID, 'units', 'meters')) 1526 call nf90_err(nf90_put_att(ncid, levID, 'positive', 'up')) 1527 call nf90_err(nf90_put_att(ncid, levID, 'standard_name', 'height')) 1528 1529 if (.not.allocated(coord)) allocate(coord(numxgrid)) 1530 do i = 1,numxgrid 1531 coord(i) = outlon0 + (i-0.5)*dxout 1532 enddo 1533 call nf90_err(nf90_put_var(ncid, lonID, coord(1:numxgrid))) 1534 deallocate(coord) 1535 if (.not.allocated(coord)) allocate(coord(numygrid)) 1536 do i = 1,numygrid 1537 coord(i) = outlat0 + (i-0.5)*dyout 1538 enddo 1539 call nf90_err(nf90_put_var(ncid, latID, coord(1:numygrid))) 1540 deallocate(coord) 1541 call nf90_err(nf90_put_var(ncid, levID, outheight(1:numzgrid))) 1542 1543 ! write time, one field per time - different to the others! 1544 call nf90_err(nf90_put_var( ncid, timeID, itime, (/ 1 /))) 1545 1546 dIDs = (/ londimid, latdimid, levdimid, timedimid, pointspecdimid, nageclassdimid /) 1547 1548 do ks=1,nspec 1549 do kp=1,maxpointspec_act 1550 do nage=1,nageclass 1551 1552 write(anspec,'(i3.3)') ks 1553 1554 ! East Flux 1555 call nf90_err(nf90_def_var(ncid,'flux_east_'//anspec, nf90_float, dIDs, & 1556 fluxID)) 1557 1558 do jy=0,numygrid-1 1559 do ix=0,numxgrid-1 1560 do kz=1, numzgrid 1561 grid(ix,jy,kz)=flux(1,ix,jy,kz,ks,kp,nage) 1562 end do 1563 end do 1564 end do 1565 1566 call nf90_err(nf90_put_var(ncid,fluxid,1.e12*grid(0:numxgrid-1,0:numygrid-1,1:numzgrid)& 1567 /areaeast(0:numxgrid-1,0:numygrid-1,1:numzgrid)/loutstep,& 1568 (/ 1,1,1,tpointer,kp,nage /), (/ numxgrid,numygrid,numzgrid,1,1,1 /) )) 1569 1570 end do 1571 end do 1572 end do 1573 1574 ! Close netCDF file 1575 call nf90_err(nf90_close(ncid)) 1576 1577 end subroutine fluxoutput_netcdf 1578 1426 1579 end module netcdf_output_mod 1427 1580 -
src/timemanager.f90
r3481cc1 rffbe224 82 82 #ifdef USE_NCF 83 83 use netcdf_output_mod, only: concoutput_netcdf,concoutput_nest_netcdf,& 84 &concoutput_surf_netcdf,concoutput_surf_nest_netcdf 84 &concoutput_surf_netcdf,concoutput_surf_nest_netcdf,fluxoutput_netcdf 85 85 #endif 86 86 … … 429 429 endif 430 430 if ((iout.eq.4).or.(iout.eq.5)) call plumetraj(itime) 431 if (iflux.eq.1) call fluxoutput(itime) 431 if ((iflux.eq.1).and.(lnetcdfout.eq.0)) call fluxoutput(itime) 432 if ((iflux.eq.1).and.(lnetcdfout.eq.1)) call fluxoutput_netcdf(itime) 432 433 write(*,45) itime,numpart,gridtotalunc,wetgridtotalunc,drygridtotalunc 433 434 -
src/timemanager_mpi.f90
r3481cc1 rffbe224 536 536 endif 537 537 if ((iout.eq.4).or.(iout.eq.5)) call plumetraj(itime) 538 if (iflux.eq.1) call fluxoutput(itime) 538 if ((iflux.eq.1).and.(lnetcdfout.eq.0)) call fluxoutput(itime) 539 if ((iflux.eq.1).and.(lnetcdfout.eq.1)) call fluxoutput_netcdf(itime) 539 540 if (mp_measure_time) call mpif_mtime('iotime',1) 540 541
Note: See TracChangeset
for help on using the changeset viewer.