[MITgcm-devel] exf_interp and grib data
Jean-Michel Campin
jmc at ocean.mit.edu
Tue Apr 14 22:30:57 EDT 2015
Hi Martin,
I don't have time to spent on this; I can see two options:
1) fixing the current exf_interp.F for the case lat_inc < 0 ;
would need to fix the special treatment of N & S poles.
but could leave exf_interp_uv.F unchanged with a stop there
if one try to use it with lat_inc < 0.
The thing that I remember when adding some of the pole treatment there
in exf_interp_uv.F is that there was several cases to check (depending on
resolution ratio) and Oliver kindly spent some time to test these; fixing
this routine for the case lat_inc < 0 would involve going back to these
test-case to check.
2) go back to original exf_interp.F (would be safer to add a stop
somewhere if lat_inc < 0 i; may be also if lon_inc < 0)
Cheers,
Jean-Michel
On Wed, Apr 08, 2015 at 09:53:10AM +0200, Martin Losch wrote:
> Hi Jean-Michel,
>
> I can reproduce your results exactly.
>
> Do we want to spend more time on this? I can simply remove the bad code and revert to the previous version. I am no longer sure how useful this ???feature??? is, given that in preparing the forcing data, one has to do some preprocessing anyway, so that flipping the y-axis of the forcing data is probable the last and least of all issues.
>
> Martin
>
> > On 08 Apr 2015, at 06:29, Jean-Michel Campin <jmc at ocean.mit.edu> wrote:
> >
> > Hi Martin,
> >
> > Did my homework, but would be good if you can check (independently).
> > It's mainly to support my second point (not sure about the logic).
> >
> > A simple set-up Nx=1, Ny=8, closed at N & S end, and testing exf-interp
> > with just one field (qnet);
> > I have ygOrigin=6., and dySpacing=2., which gives me yC =[7:2:21]
> >
> > 1rst case (lat_inc >0 )
> > since I specify:
> > hflux_lat0 = 6.,
> > hflux_lat_inc = 4., 8., 4.,
> > I assume input data is at lat: 6. 10. 18. 22.
> > and the qnet values are: +2, +2, -1, -1,
> >
> > 2nd case (lat_inc < 0)
> > hflux_lat0 = 22.,
> > hflux_lat_inc = -4., -8., -4.,
> > I assume input data is at lat: 22. 18. 10. 6.
> > and I did flip the exf input file.
> >
> > I have attached a plot, after 1 iter, of diagnostics EXFqnet
> > in the 2 experiments (and in black the original input data):
> >> fig_exf_iterp.jpg.gz
> >
> > and also the data.exf corresponding to these 2 case (data.exfi.dy+ & data.exf.dy-)
> > with my two binary input files (qnet.bin & qinv.bin).
> >
> > Cheers,
> > Jean-Michel
> >
> >
> > On Tue, Apr 07, 2015 at 09:23:50AM -0400, Jean-Michel Campin wrote:
> >> Hi Martin,
> >>
> >> I don't have time right now, but will take a look at this sometime this week.
> >>
> >> Cheers,
> >> Jean-Michel
> >>
> >> On Tue, Apr 07, 2015 at 11:53:39AM +0200, Martin Losch wrote:
> >>> Hi Jean-Michel,
> >>>
> >>> what should be do about this code?
> >>>
> >>> I think your concern (2) is OK, because the roles of s_ind and w_ind are exchanged here and s_ind is the index in the north (not the south anymore).
> >>>
> >>> The pole business looks very complicated to me, especially in conjunction with the vector interpolation in exf_interp_uv. I am wondering if the effort in coding this is maybe far larger than the effort in swapping the input fields from the atmopheric convention to ours? What do you think?
> >>> I could simply revert to the previous version and remove my addtions.
> >>>
> >>> Martin
> >>>
> >>>> On 24 Mar 2015, at 19:18, Jean-Michel Campin <jmc at ocean.mit.edu> wrote:
> >>>>
> >>>> Hi Martin,
> >>>> Enjoy your vacation.
> >>>> Cheers,
> >>>> Jean-Michel
> >>>>
> >>>> On Tue, Mar 24, 2015 at 07:00:32PM +0100, Martin Losch wrote:
> >>>>> Hi Jean-Michel,
> >>>>>
> >>>>> I guess you are right and my code does not really work properly, at least I totally forgot the UV part.
> >>>>>
> >>>>> 2. ind_s and w switch their meaning for latinc <0, so that part should be OK, but not sure right now.
> >>>>>
> >>>>> I am on the run to catch the vacation train. Lets sort that out when I return.
> >>>>>
> >>>>> Martin
> >>>>>
> >>>>>> On 24.03.2015, at 16:27, Jean-Michel Campin <jmc at ocean.mit.edu> wrote:
> >>>>>>
> >>>>>> Hi Martin,
> >>>>>>
> >>>>>> I've just took a quick look at the changes you made in exf_interp.F
> >>>>>> and I have few questions:
> >>>>>> 1) I am not sure about the poles treatment (from line 139 to 175)
> >>>>>> when lat_inc is < 0.
> >>>>>> 2) and I am not sure I am fully getting the logic:
> >>>>>> the model grid point lat is found to be between y_in(s_ind)
> >>>>>> and y_in(s_ind-1) (at least it's what we check on line 412, right ?)
> >>>>>> But in S/R EXF_INTERPOLATE, we use s_ind and s_ind+1 values,
> >>>>>> so I am confused.
> >>>>>> Or is it that s_ind is no longer the index of the data point
> >>>>>> that is just south of the model point ?
> >>>>>> The experiment global_with_exf is great but it does not really
> >>>>>> use the interpolation coeff since the data grid and model grid
> >>>>>> are the same.
> >>>>>> 3) seem that S/R EXF_INTERP_UV should also have been modified (or
> >>>>>> put some check and stop there) in case lat_inc < 0.
> >>>>>>
> >>>>>> Cheers,
> >>>>>> Jean-Michel
> >>>>>>
> >>>>>>> On Fri, Mar 13, 2015 at 03:00:16PM +0100, Martin Losch wrote:
> >>>>>>> Yes, all in global_with_exf/input
> >>>>>>>
> >>>>>>> Martin
> >>>>>>>
> >>>>>>>> On 13 Mar 2015, at 14:45, Jean-Michel Campin <jmc at ocean.mit.edu> wrote:
> >>>>>>>>
> >>>>>>>> Hi Martin,
> >>>>>>>>
> >>>>>>>> You can go for it. The new wind-stress file could stay in
> >>>>>>>> global_with_exf/input
> >>>>>>>> does not need to be in tutorial_global_oce_latlon/input, right ?
> >>>>>>>>
> >>>>>>>> Cheers,
> >>>>>>>> Jean-Michel
> >>>>>>>>
> >>>>>>>>> On Thu, Mar 12, 2015 at 05:02:47PM +0100, Martin Losch wrote:
> >>>>>>>>> I think I found a way to do this. Since it is not so difficult after all, I???d like to check in the code.
> >>>>>>>>>
> >>>>>>>>> We could test this with some (maybe just one?) of the files in global_with_exf, which is the only experiment that tests the interpolation code anyway.
> >>>>>>>>> I can get perfect agreement with ???trenberth_taux.bin" flipped (I call this now trenberth_taux.flipped???), and
> >>>>>>>>> ustress_lat0 = 78.,
> >>>>>>>>> ustress_lat_inc = 39*-4.,
> >>>>>>>>>
> >>>>>>>>> The only problem with this test is, that it will increase the size of the MITgcm download by the size of trenberth_taux.bin (well, 169kB). What do you think?
> >>>>>>>>>
> >>>>>>>>> Martin
> >>>>>>>>>
> >>>>>>>>>
> >>>>>>>>>
> >>>>>>>>>> On 12 Mar 2015, at 15:26, Martin Losch <Martin.Losch at awi.de> wrote:
> >>>>>>>>>>
> >>>>>>>>>> Hi Jean-Michel,
> >>>>>>>>>>
> >>>>>>>>>> I tried to add code for the case of y_in(0)>y_in(nyIn+1) (or lat_inc<0), but wasn???t successful.
> >>>>>>>>>> I guess I need to understand this code first, but before I seriously try, do you think it???s worth it?
> >>>>>>>>>>
> >>>>>>>>>> Is this the only place where this assumption is made? Or maybe there???s something exf_interpolate.F as well?
> >>>>>>>>>>
> >>>>>>>>>> Martin
> >>>>>>>>>>
> >>>>>>>>>>> On 12 Mar 2015, at 14:15, Jean-Michel Campin <jmc at ocean.mit.edu> wrote:
> >>>>>>>>>>>
> >>>>>>>>>>> Hi Martin,
> >>>>>>>>>>>
> >>>>>>>>>>> I did not check carefully, but I vagely remember that exf-interp assumes
> >>>>>>>>>>> lat_inc > 0 (& lon_inc > 0). For instance, exf_interp.F, lines 350-357:
> >>>>>>>>>>> IF ( w_ind(i,j).GT.s_ind(i,j)+1 ) THEN
> >>>>>>>>>>> k = NINT( (s_ind(i,j)+w_ind(i,j))*0.5 )
> >>>>>>>>>>> IF ( yG(i,j,bi,bj) .LT. y_in(k) ) THEN
> >>>>>>>>>>> w_ind(i,j) = k
> >>>>>>>>>>> ELSE
> >>>>>>>>>>> s_ind(i,j) = k
> >>>>>>>>>>> ENDIF
> >>>>>>>>>>> ENDIF
> >>>>>>>>>>> In theory, should be possible to (re-)write the code for the case
> >>>>>>>>>>> of negative increment but, except if I am wrong, this is not presently the case.
> >>>>>>>>>>>
> >>>>>>>>>>> Cheers,
> >>>>>>>>>>> Jean-Michel
> >>>>>>>>>>>
> >>>>>>>>>>>> On Thu, Mar 12, 2015 at 11:18:41AM +0100, Martin Losch wrote:
> >>>>>>>>>>>> Hi there
> >>>>>>>>>>>>
> >>>>>>>>>>>> I am trying to convert (on the fly) atmospheric data in grib format into ieee-be so that I can use exf (exf_interp) to read it.
> >>>>>>>>>>>> The grb data seems to put the index (i,j) = (1,1) at the top left (NW) corner, whereas we expect (i,j)=(1,1) to be at the bottom right (SW) corner of the gridded data field. For technical reasons, it would be simpler to not have to flip the forcing fields (I would like to use wgrib <http://www.cpc.ncep.noaa.gov/products/wesley/wgrib.html>, which can produce ieee-be fields from grib data, but doesn???t seem to be able to do the N to S-flip), so that my data.exf entry would look like this:
> >>>>>>>>>>>>
> >>>>>>>>>>>> ${var}_lat0 = 78,
> >>>>>>>>>>>> ${var}_lat_inc = 39*-4.,
> >>>>>>>>>>>>
> >>>>>>>>>>>> instead of
> >>>>>>>>>>>> ${var}_lat0 = -78,
> >>>>>>>>>>>> ${var}_lat_inc = 39*4.,
> >>>>>>>>>>>>
> >>>>>>>>>>>> Appart from the fact, that the routine exf_interp checks if the forcing grid ???encompasses??? the model grid (which it does technically), and this fails (because the check also expects lat0 to be the southern-most latitude), can exf_interp handle this situation properly? (a simple test with global_with_exf suggests, that it does not, but why?)
> >>>>>>>>>>>>
> >>>>>>>>>>>> Martin
> >>>>>>>>>>>>
> >>>>>>>>>>>>
> >>>>>>>>>>>>
> >>>>>>>>>>>> _______________________________________________
> >>>>>>>>>>>> MITgcm-devel mailing list
> >>>>>>>>>>>> MITgcm-devel at mitgcm.org
> >>>>>>>>>>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
> >>>>>>>>>>>
> >>>>>>>>>>> _______________________________________________
> >>>>>>>>>>> MITgcm-devel mailing list
> >>>>>>>>>>> MITgcm-devel at mitgcm.org
> >>>>>>>>>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
> >>>>>>>>>>
> >>>>>>>>>>
> >>>>>>>>>> _______________________________________________
> >>>>>>>>>> MITgcm-devel mailing list
> >>>>>>>>>> MITgcm-devel at mitgcm.org
> >>>>>>>>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
> >>>>>>>>>
> >>>>>>>>>
> >>>>>>>>> _______________________________________________
> >>>>>>>>> MITgcm-devel mailing list
> >>>>>>>>> MITgcm-devel at mitgcm.org
> >>>>>>>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
> >>>>>>>>
> >>>>>>>> _______________________________________________
> >>>>>>>> MITgcm-devel mailing list
> >>>>>>>> MITgcm-devel at mitgcm.org
> >>>>>>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
> >>>>>>>
> >>>>>>>
> >>>>>>> _______________________________________________
> >>>>>>> MITgcm-devel mailing list
> >>>>>>> MITgcm-devel at mitgcm.org
> >>>>>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
> >>>>>>
> >>>>>> _______________________________________________
> >>>>>> MITgcm-devel mailing list
> >>>>>> MITgcm-devel at mitgcm.org
> >>>>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
> >>>>>
> >>>>> _______________________________________________
> >>>>> MITgcm-devel mailing list
> >>>>> MITgcm-devel at mitgcm.org
> >>>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
> >>>>
> >>>> _______________________________________________
> >>>> MITgcm-devel mailing list
> >>>> MITgcm-devel at mitgcm.org
> >>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
> >>>
> >>>
> >>> _______________________________________________
> >>> MITgcm-devel mailing list
> >>> MITgcm-devel at mitgcm.org
> >>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
> > <fig_exf_iterp.jpg.gz><data.exf.dy+><data.exf.dy-><qnet.bin><qinv.bin>_______________________________________________
> > MITgcm-devel mailing list
> > MITgcm-devel at mitgcm.org
> > http://mitgcm.org/mailman/listinfo/mitgcm-devel
>
>
> _______________________________________________
> MITgcm-devel mailing list
> MITgcm-devel at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-devel
More information about the MITgcm-devel
mailing list