Odd values from SG smoother

[Dave Offiler]

Hi Dave,

when running your t_roppthin.sh script, I get differences between the unthinned and thinned profile. For REF and I think also BA, there is a large deviation at the lowest altitude for the SG cases. For the longs, there is one value different to the other ones at the lowest altitude. And sorry, it is the longs, not the lats.

I noticed something in my tests where I suddenly found a very low geop value in the thinned data, although the unthinned data wasn't covering that range. So I assume it is the SG processing.

I haven't done any further testing, but it should be easy for you to at least see the problem.

Cheers,

Axel

[Dave Offiler]

The stand-alone SG smoother appears to have problems with both first & last values in the smoothed array, in that the end points do not have valid neighbours to use in the averaging kernel.

As a temporary fix, I've restricted the returned smoothed array to points 2 to nLev-1. This removes the problem values at each end of the array (though only the bottom of the profile in the full application seemed to be in error.) The thinning test does not now show the problem values at the bottom of the profiles - ifort check only.

This code change committed as [1096]

To be properly fixed:

1) the 'good' zone should be nl+1 to nlev-nr in case the SG parameters are ever changed (both nl [no. pts to left] & nr [no. pts. to right] are currently set to 1)

2) the multiplication in the BLAS routine dgbmv may be trying to use array elements 1-nl and nLev+nr (ie going out of bounds) when calculating a value for the true end points. So using whatever random value happens to be in memory at that point is getting into the 1st & nLev-th value. To avoid possible problems with other compilers (and even if ifort doen't crash, this isn't A Good Thing), some way of preventing out-of-bounds referencing should be sought.

[Dave Offiler]

Sean Healy rang up at around 11:40 this morning, reporting a large global bias in CHAMP/GRACE-A at ~40km vs ECMWF model BAs.

This may well be a manifestaton of the SG problem at the top end of the original GFZ profiles (not at the top end of the thinned levels). I have now done a 'make install' and copied the executable to my own bin directory so this fixed version will be used from (just before) mid-day for the GTS encoding. Sean will report back tomorrow on whether the bias has been removed or not after this time.

[Dave Offiler]

Took apart the BLAS dgbmv code and re-coded the (relatively trival) matrix multiplication in-line to see what was going on. The code doen't go out of array bound (proven with "ifort -check bounds") but the inner loop does get clamped at the ends, which is likely cause for the 1st & last elements of the array to be not correctly calculated. The returned portion of the smoothed array now uses the nl & nr values rather than hard-coded '1' values.

Rechecked using: 1) a stand-alone call to ropp_io_thin_sg() - creating a dummy N profile with 10% random noise and outputting before-and-after profiles and plotting both in PV-WAVE 2) the full ropp2ropp test in script t_roppthin.sh Smoothed profiles now compare as expected with original profile.

Committed as part of [1102].

Leaving open pending Sean's response.

[Dave Offiler]

Sean has not explicitly checked the bias (ECMWF are not using GFZ data operationally), but doesn't see any problems. Closing this ticket as fixed.