Search for 'longitude' and 'latitude' variables when reading .nc DEM files

[chlnddev]

Search for 'longitude' and 'latitude' variables when reading .nc DEM files

[krober10nd]

Hey Ravi, Thanks for your help and support. I will take care of this later today. We will need to check a cell array of known nonmenclatures that describe coordinates in DEMs using regexp.

…

On Tue, Jun 4, 2019, 12:13 PM Javi RG ***@***.***> wrote: (Sorry for being so annoying - I only want to improve this project) When the program reads the netCDF DEM file, only checks for the existence of 'lon', 'lat', 'x' or 'y' variables. In some datasets (including the SRTM 15 Plus), those variables are called 'longitude' and 'latitude'. geodata.m at line 161 must try to find those nomenclatures before throwing an error. — You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub <#8?email_source=notifications&email_token=AEOBZ7FZ6R2DQ4Z47OHVC6TPY2IDJA5CNFSM4HS6OU5KYY3PNVWWK3TUL52HS4DFUVEXG43VMWVGG33NNVSW45C7NFSM4GXSJO5A>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AEOBZ7HV2CYDYQ6XYDVAQZDPY2IDJANCNFSM4HS6OU5A> .

Most of them follow the CF conventions, so the array won't be huge: http://cfconventions.org/cf-conventions/cf-conventions.html#coordinate-types

[krober10nd]

I've updated the Projection branch to support the three basic x- and y- coordinate variables names in DEM data that is CF compliant (see the arrays wkv_x and wkv_y).

Note also that I've refactored the geodata class so please let me know if things have stopped working the same way they were!

Keith,

Two things to comment.

1. Longitude and Latitude words on lines 49 & 50 of geodata.m must be lowercase according to CF convention, or include the lowercase ones in the array as well as the title-case words.

2. This part of geodata.m is not a good piece of code (line 400 and beyond)
   Once x & y values are successfully read, you should not try to read again those values from the DEM file.

for i = 1 : length(obj.wkv_x)
    try
        x = double(ncread(obj.demfile,obj.wkv_x{i}));
        y = double(ncread(obj.demfile,obj.wkv_y{i}));
    catch
        if i == length(obj.wkv_x)
           error('Could not locate x/y coordinate in DEM') ; 
        end
    end
end

I'm not so keen on Matlab because I work in C++, but something like this should solve the issue:

if ~isempty(x) && ~isempty(y)
    break % Get out of the loop
end 

[krober10nd]

Hello Javi,

The try and catch structure is necessary in MATLAB to handle the error produced by the NetCDF wrapper.

The x and y coordinates are read twice to handle a variety of use cases (i..e, no bounding box provided, no DEM, etc.) to exit quickly. The code below 400 is to deal with the situation of a regional domain straddling the 180 meridian, which turns out to be a non-trivial problem for automatic mesh generation.

If you have code changes you'd like to incorporate please submit a pull request. Thanks.

All the code under the if ~isempty(fname) only affects when using a DEM file. The x and y should be read as many times as needed until we match the variable in the netCDF.

Right now, your code is reading the correct values and then it overwrites them. That is not the expected behavior.

[WPringle]

I reopened issue because I also believe the code is as was faulty. I just made a commit that should hopefully fix this issue. @JaviRG can you please check it works for you.

[krober10nd]

What precisely was faulty?

[WPringle]

As noted by @JaviRG

Right now, your code is reading the correct values and then it [attempts to] overwrites them. That is not the expected behavior.

The break out of the loop once x and y are correctly read is better to avoid unnecessary read attempts while relying on the try and catch structure.
Perhaps faulty is not correct because it worked but my commit was aimed as at trying to avoid a clumsy try and catch and instead to search the information that is available in netcdf for reading.
Note that actual bugs were fixed too since 'fname variable' was used as the netcdfile name in ParseDEM but obj.demfile was then used for reading x and y variables, along with another ubiquitous == 0 versus isempty() catch issue.

[krober10nd]

Fair enough. In the future, the person should submit a pull request and we inspect the code. I do not want to commit un tested code into the branch that's supposed to always work and then rely on users to verify it's working as expected.

[WPringle]

That's true. After testing our common examples I believe that the most recent commit should work.
As mentioned by @krober10nd we encourage a pull request for suggested code changes to accompany an issue if a plausible fix is known to the OP of the issue.

Thanks @WPringle for reopening the issue and creating a new function to handle this. Now it's working properly.

In the future, the person should submit a pull request and we inspect the code.

I would be glad to contribute writing code, but as I mentioned before, I work with C++ and don't know Matlab. I think (at least) is better to open an issue and tell the authors about weaknesses in the code.