martagilbardaji
New member
Hello all,
I am testing meshes that have a central region with spacings of around 1km, and my problem is that, for some of the cells in the central fine-resolution region, the terrain ('ter') value given by init_atmosphere when computing the static.nc file is NaN. However, the static.nc file is generated without throwing errors, and I only realised this was happening because using that static.nc file to generate the initial conditions produced a segmentation fault due to some zgrig(k, iCell) being NaNs.
The issue
In this picture, the terrain ('ter') variable of the static.nc file is plotted (zooming over the central fine-resolution region). There are 7 cells that have invalid (NaN) terrain values, colored dark grey.
I attach a zip folder (terrainNans.zip) with the input regional grid.nc file; the namelist.init_atmosphere and streams.init_atmsphere files I used; and the output static.nc file that has been generated. The WPS_GEOG files are used were downloaded from the offical webpage, and from a message here in the forum that pointed to some extra files (landuse_30s and modis_landuse_20class_30s I think).
My theory
I have checked the source code and it seems to me that the HGT interpolation in mpas_init_atm_static.F does not check whether there are cells that can end up without terrain information. As I understand, the interpolation procedure consists of assigning to each MPAS cell the mean terrain of the topo_gmted2010_30s grid-points contained in that MPAS cell. However it does not take into account that it is possible that some MPAS cells, if they have fine grid-spacings, do not contain any of the topo_gmted2010_30s grid-points. Those cells end up with NaN values and no error is detected. I think that there should be a at least test for that situation, and it would be nice to implement a solution/algorithm that avoids this issue.
For now, in my source code, I added this snippet :
(...)
do iCell = 1,nCells
!MGB Test that no cells are empty of terrain data
if (nhs(iCell) == 0) then
call mpas_log_write('--- Cell $i was not filled by any topo', intArgs=(/iCell/))
call mpas_log_write('*****************************************************************', messageType=MPAS_LOG_ERR)
call mpas_log_write('Invalid topography dataset: at least one cell has NaN terrain ', messageType=MPAS_LOG_ERR)
call mpas_log_write(' This happens because the grid-spacing is finer than in the reference topo dataset.', messageType=MPAS_LOG_ERR)
call mpas_log_write('*****************************************************************', messageType=MPAS_LOG_ERR)
call mpas_log_write('Use finer reference topo or modify the code.', messageType=MPAS_LOG_CRIT)
end if
ter(iCell) = ter(iCell) / real(nhs(iCell))
end do
deallocate(rarray)
deallocate(nhs)
call mpas_log_write('--- end interpolate TER')
(...)
In fact, this also happens for other variables like the vegetation fraction (greenfrac) variable and soil category (isltyp), and the incorrect value is 0 instead of NaN, which is even worse because the next steps of the process will probably not fail but will be using very wrong data. So, I think tests should be added after all the interpolations that can potentially have this issue.
Solutions?
A solution I guess it would be to create a higher-resolution topography/landuse... option for the WPS input files, but I do believe something more elegant can be done that does not invlove using up much more storage space.
What do you suggest? Has someone found a work-around this issue? I believe some studies have used MPAS with resolutions around 1km, so I assume they managed to solve this, but I could not find anything here in the forum.
Any help will be appreciated,
Marta
I am testing meshes that have a central region with spacings of around 1km, and my problem is that, for some of the cells in the central fine-resolution region, the terrain ('ter') value given by init_atmosphere when computing the static.nc file is NaN. However, the static.nc file is generated without throwing errors, and I only realised this was happening because using that static.nc file to generate the initial conditions produced a segmentation fault due to some zgrig(k, iCell) being NaNs.
The issue
In this picture, the terrain ('ter') variable of the static.nc file is plotted (zooming over the central fine-resolution region). There are 7 cells that have invalid (NaN) terrain values, colored dark grey.
I attach a zip folder (terrainNans.zip) with the input regional grid.nc file; the namelist.init_atmosphere and streams.init_atmsphere files I used; and the output static.nc file that has been generated. The WPS_GEOG files are used were downloaded from the offical webpage, and from a message here in the forum that pointed to some extra files (landuse_30s and modis_landuse_20class_30s I think).
My theory
I have checked the source code and it seems to me that the HGT interpolation in mpas_init_atm_static.F does not check whether there are cells that can end up without terrain information. As I understand, the interpolation procedure consists of assigning to each MPAS cell the mean terrain of the topo_gmted2010_30s grid-points contained in that MPAS cell. However it does not take into account that it is possible that some MPAS cells, if they have fine grid-spacings, do not contain any of the topo_gmted2010_30s grid-points. Those cells end up with NaN values and no error is detected. I think that there should be a at least test for that situation, and it would be nice to implement a solution/algorithm that avoids this issue.
For now, in my source code, I added this snippet :
(...)
do iCell = 1,nCells
!MGB Test that no cells are empty of terrain data
if (nhs(iCell) == 0) then
call mpas_log_write('--- Cell $i was not filled by any topo', intArgs=(/iCell/))
call mpas_log_write('*****************************************************************', messageType=MPAS_LOG_ERR)
call mpas_log_write('Invalid topography dataset: at least one cell has NaN terrain ', messageType=MPAS_LOG_ERR)
call mpas_log_write(' This happens because the grid-spacing is finer than in the reference topo dataset.', messageType=MPAS_LOG_ERR)
call mpas_log_write('*****************************************************************', messageType=MPAS_LOG_ERR)
call mpas_log_write('Use finer reference topo or modify the code.', messageType=MPAS_LOG_CRIT)
end if
ter(iCell) = ter(iCell) / real(nhs(iCell))
end do
deallocate(rarray)
deallocate(nhs)
call mpas_log_write('--- end interpolate TER')
(...)
In fact, this also happens for other variables like the vegetation fraction (greenfrac) variable and soil category (isltyp), and the incorrect value is 0 instead of NaN, which is even worse because the next steps of the process will probably not fail but will be using very wrong data. So, I think tests should be added after all the interpolations that can potentially have this issue.
Solutions?
A solution I guess it would be to create a higher-resolution topography/landuse... option for the WPS input files, but I do believe something more elegant can be done that does not invlove using up much more storage space.
What do you suggest? Has someone found a work-around this issue? I believe some studies have used MPAS with resolutions around 1km, so I assume they managed to solve this, but I could not find anything here in the forum.
Any help will be appreciated,
Marta