to_latlon does not rotate velocities

[dhruvbalwada]

This is related to #204 .

One might think that faces_dataset_to_latlon would rotate the velocities. But it does not seem to.

This an example of what the meridional velocity looks like in the North Atlantic (generated in https://github.com/dhruvbalwada/LLC_analysis/blob/master/LLC_grid_visualize.ipynb - see near bottom). The Gulf Stream shows up as a strong meridional band, rather than a zonal current.

Has anyone managed to fix this?

Also, is there a simple way to make the plot so that there is no banding where the faces meet?

[rabernat]

That looks like a bug. It's not supposed to do that.

[rabernat]

The rotation happens here

https://github.com/MITgcm/xmitgcm/blob/f5ee774f6ee3788d48f32a9bd1fbe02ebea9449b/xmitgcm/llcreader/llcmodel.py#L231-L234

[rabernat]

Just as a verification that this rotation is happening correctly, here is a quick example with ECCO (LLC90)

from intake import open_catalog
from xmitgcm import llcreader

cat = open_catalog("https://raw.githubusercontent.com/pangeo-data/pangeo-datastore/master/intake-catalogs/ocean.yaml")
ds  = cat["ECCOv4r3"].to_dask()
ds_ll = llcreader.llcmodel.faces_dataset_to_latlon(ds)
ds_ll.UVELMASS[0, 0].load().plot(figsize=(20, 10))

There is a lot of other stuff going on with your example (xgcm, custom plotting, etc.), so I would recommend trying to simplify things down to get to the source of the problem. When debugging llcreader, I highly recommend working with an llc90 dataset, as it is much faster to compute everything. All of the issues are (or at least shoudl be) the same.

One thought: the weird artifacts of your plot are indicative of what happens when the x and y variables for pcolormesh have missing data. You are plotting like this

plt.pcolormesh(V_sel_coarse.XC, V_sel_coarse.YC, V_sel_coarse)

Are you sure there are no NaNs in V_sel_coarse.XC and V_sel_coarse.YC? It is unfortunate that MITgcm outputs this grid data with missing values, but that's beyond our control and we have to work around it somehow.

[rabernat]

@dhruvbalwada - any update on this? Have you checked my suggestion?

Are you sure there are no NaNs in V_sel_coarse.XC and V_sel_coarse.YC?

[rabernat]

@dhruvbalwada I'd like to resolve this. Did you get a chance to check my suggestion?

[dhruvbalwada]

I will take a look at this this week, been a bit busy.

[dhruvbalwada]

I tried your above suggestion: https://gist.github.com/dhruvbalwada/a2800d9fa0ee1083d9b064a7f93105e9 It seems that for the 4320 the automatic rotation is not taking place for some reason, unless the way I am loading the data (ln:13) for the 4320 is not a compatible way to do things. Notice how for the low res model the coordinates are automatically loaded, but i needed to do it manually for 4320.

[dhruvbalwada]

I managed to fix the problem with a hack. The problem it seems is that the SSU and SSV velocities in the catalogue do not have the appropriate attributes, and so the faces_dataset_to_latlon function is unable to figure out that the U and V are "mate" variables.

Adding the appropriate mate attribute fixed the issue:

u = cat_4320.LLC4320_SSU.to_dask()
v = cat_4320.LLC4320_SSV.to_dask()
coords  = cat_4320.LLC4320_grid.to_dask()

u.U.attrs = {'long_name': 'Zonal Velocity (m/s)', 
            'mate': 'V', 'units':'m/s'}
v.V.attrs = {'long_name': 'Meridional Velocity (m/s)', 
            'mate': 'U', 'units':'m/s'}

ds_4320 = xr.merge([u, v, coords])

Here is a plot of the zonal velocity component (with slight coarsening) :

Any suggestion on how to resolve this issue for everyone? Can the attributes be added to the data on the intake catalogue?

[csy0101]

我设法通过黑客解决了这个问题。 看来问题在于目录中的SSU和SSV速度没有适当的属性，因此该faces_dataset_to_latlon函数无法确定U和V是“匹配”变量。

添加适当的mate属性可解决此问题：

u = cat_4320.LLC4320_SSU.to_dask()
v = cat_4320.LLC4320_SSV.to_dask()
coords  = cat_4320.LLC4320_grid.to_dask()

u.U.attrs = {'long_name': 'Zonal Velocity (m/s)', 
            'mate': 'V', 'units':'m/s'}
v.V.attrs = {'long_name': 'Meridional Velocity (m/s)', 
            'mate': 'U', 'units':'m/s'}

ds_4320 = xr.merge([u, v, coords])

这是纬向速度分量的图（略有粗化）：

关于如何为每个人解决此问题的任何建议？可以将属性添加到摄入量目录中的数据中吗？

My original work plan was to use the u, v of the llc4320 to compare the sea surface velocity retrieved by satellites, but now I have a question, are these two types of u, v the same property?

[rabernat]

Any suggestion on how to resolve this issue for everyone? Can the attributes be added to the data on the intake catalogue?

This is great Dhruv! I think we just need to add a bit of metadata to the zarr stores for LLC4320. Good job digging to the bottom of this! 😄

[dhruvbalwada]

What do we need to do close this? Should we tag Charles to add the metadata?

[rabernat]

There are some subtleties to deal with here (consolidated metadata, etc.), so I want to do it myself. But I have some urgent deadlines this week, and it is not a high priority right now. The workaround documented in this issue should suffice for anyone who needs to do this within the next week or so.

Once I update the metadata, I'll close the issue.