sigma of LSC2 class

[feiye-vims]

@felicio93 @josephzhang8

Hey Felicio,

We noticed some issues with the new LSC2 class when running our existing scripts for generating nudging files.

The issue is that some scripts rely on the "sigma" attribute of an instance of the LSC2 class. I implemented a temporary fix. You can review it and make additional adjustment to the class as you see fit.

I see that you used a few parameters in "def init(self, hsm, nv, h_c, theta_b, theta_f)" and then calculate the sigma layers ("_snd"), which is a valid procedure. However, we also want to be able to read an existing vgrid.in without knowing how it is made (i.e., without knowing the parameter list of init ), and access something like "my_vgrid.sigma" later.

Thanks!

[felicio93]

Sounds good, let me take a look and I will get back to you soon.

[felicio93]

I am not sure what has changed but to use pyschism I had to:
Create an environment with python 3.9 (3.10 did not work for me) and install the following libraries:

conda create -n pyschism python=3.9
pip install pyschism
conda install -c conda-forge cfunits
pip install fiona

[feiye-vims]

I'm not sure either. I have pyschism with python 3.10 on SciClone, but I did the installation months ago.
I'll install it on Hercules inside the "stofs" environment (also based on 3.10) next week for the 2D setup, let me see if I come across the same problem.

[felicio93]

Hi Fei, I reviewed the code. Tell me if I understand what you mean.

Now one can open an vgrid.in file with:

v_grid = LSC2.open(vgrd_path)

and return the sigma with: v_grid.sigma. One can also return nv with v_grid.nv, which is being defined by from_sigma.

Are you saying you also need/want it to also return hsm, h_c, theta_b, theta_f? or we are good with returning the sigma values only? Please let me know if that is the case, I see that the SZ class already does something similar, maybe I can borrow things from it. After giving more thought I am not sure if it is possible to go from a lsc2 vgrid.in to the other variables (hsm, h_c, theta_b, theta_f). pyschism does that for the SZ-based vgrid.in by parsing the vgrid file once all variables are written within it. That is not the case for the LSC2-based vgrid.in, right?

by the way, I don't think nv = sigma.shape[1] is correct. nv is an array that contains the number of layers per depth threshold (hsm). I think sigma.shape[1] is nvrt not nv.

[feiye-vims]

Yes, I need LSC2.open to read an existing vgrid.in. I think this function has been there for a while. It didn't work yesterday because it tries to return cls(sigma). When reading an existing vgrid.in, we do not need to know or return hsm, h_c, theta_b, theta_f. Like you said, they cannot be derived if we only have vgrid.in.

On the other hand, if we initialize LSC2 with def init(self, hsm, nv, h_c, theta_b, theta_f), it is perfectly fine to have these attributes exposed as v_grid.hsm, v_grid.h_c, etc. In addition, it is helpful to also have v_grid.nvrt and v_grid.sigma.

Yes, I meant nvrt = sigma.shape[1], please help correct it.

[felicio93]

Great,
But Fei, you can do that already, right?

if you open a lsc2-based vgrid.in and then call the sigma it will return what you want (at least what I think you want), for instance:

lsc2_obj_open = LSC2.open(r'PATH\vgrid_lsc2.in')
lsc2_obj_open.sigma

The value it returns is the same as _snd (but "flipped", and with -1 instead of nan).
For instance, if I create a vgrid.in, I have _snd values of:

Now, if I save this file and open it with LSC2.open:

Yes, I can correct the nvrt = sigma.shape[1], please let me know if I misunderstood what you mean/need. I will be glad to help with this!

[feiye-vims]

Yes, I committed a change yesterday to get LSC2.open working, which I referred to as "a temporary fix" in the original message. I added a factory method "from_sigma" and changed the return value from cls(sigma) to cls.from_sigma(sigma).

Since you wrote most of the content for the LSC2 class, I thought you might prefer to structure it some other way as you see fit. If you don't see anything to change, the class works fine as is (beside the nv typo).

Some additional thoughts: @josephzhang8 Please see if this makes sense.
In the future, we may also consider slightly restructuring the class to initialize an LSC2 instance with only the sigma array and nothing else. The nvrt and kbp (bottom layer index at each node) can be derived from the sigma array. To me, calculating the sigma array using "calc_lsc2_att" based on hsm, nv, theta_b, etc. is more like a factory method. There could be other methods to make an LSC2 vgrid. Then, we will need more methods like "calc_lsc2_att_from_A", "calc_lsc2_att_from_B" ... . In fact, there's nothing stopping you to specify a different vertical discretization at each mesh node by hand. In other words, the sigma array is essential to an LSC2 vgrid.in, while the other parameters are only tied to one of the methods of generating the sigma array. But let's worry about this later.

[feiye-vims]

The value it returns is the same as _snd (but "flipped", and with -1 instead of nan).
For instance, if I create a vgrid.in, I have _snd values of:

I didn't notice the flipping and I directly make a reference to sigma when initializing _snd.

Okay, let me fix this as well as the nv typo on my side.

[felicio93]

Fei feel free to change the LSC2 class anyway you see fit! I really appreciate you letting me know but feel free to change it as much as need. Also feel free to ask me to change it if you are busy with other things - I will be glad to help you!!!

[feiye-vims]

Sounds good. Sorry if the original message was confusing, just want to be polite when I change others' code.