Problems with showing geolocated data with geoshow. Error in Zdata

Hi. I have a number of netCDF files containing geolocated satellite data that I need to plot and analyze, but I’m having some difficulties doing so. For example I have a file called COTdm201101010000001230023301MA.nc containing optical cloud thickness for a specific time. This is how I try to open and read it:

netCDF_info = ncinfo('COTdm201101010000001230023301MA.nc');

lat = ncread('COTdm201101010000001230023301MA.nc','lat');

lon = ncread('COTdm201101010000001230023301MA.nc','lon');

cot = ncread('COTdm201101010000001230023301MA.nc','cot');

And at this point everything works and I have one variable “lat” with 3600 latitude values from -90° to +90°. Another one called “lon” also with 3600 values from -90° to +90° and finally a variable called “cot” with 3600x3600 values for cloud thickness.

If I understand the documentation behind the mapping toolbox correctly this should be a simple matter of using “geoshow(lat,lon,Z)”, which in this case would be “ geoshow(lat,lon,cot) ”. But when I do this I get the error:

Warning: Error creating or updating Surface
 Error in value of property  ZData
 Array is wrong shape or size

Which I don’t really understand, to me the array seems to have exactly the correct shape. What am I doing wrong here? I should also mention that there are some gaps in the cot variable that are “NaN” if I look inside the variable, could that be the problem or am I just doing something fundamentally wrong with the geoshow? Maybe it is that my data isn’t covering the entire globe and geoshow expects it to do so?

[Rob Comer]

Can you say what size(lat), size(lon), and size(cot) are?

[PetterS]

yes.

size(lat)= 3600   1
size(lon)= 3600   1
size(cot)= 3600   3600

I suspected that the problem could be that both lat and long are in the format 3600x1 and instead they should be 3600x1 and 1x3600 respectively, but I believe that I have tried to transpose them in all combinations that I can think of and it doesn’t matter.

Answer 1

[Chad Greene]

Use meshgrid.

[latgrid,longrid] = meshgrid(lat,lon); 
geoshow(latgrid,longrid,cot,'DisplayType','texturemap')

[PetterS]

Yes that worked, thank you!

[Chad Greene]

great!

Answer 2

[Kelly Luetkemeyer]

You can explicitly tell geoshow to treat your data as a surface. Try either of these commands:

    geoshow(lat,lon,cot,'DisplayType','surface')

or

   geoshow(lat,lon,cot,'DisplayType','texturemap')

[PetterS]

I tried both commands but unfortunately I’m getting exactly the same error message as before