[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: Fwd: Re: 20040322: Any idea why gribtocdl does a "floating exception" on this file?

Rich,

I tried with doubles and the end points still are not correct. Since the
true lat/lon points are mainly for human information, I was going to
create a cdl without the all the points but supply the  true lat/lon
begging and ending points. What do you think?  gribtonc will still decode
the data into netcdf.

RObb...


On Fri, 2 Apr 2004, Rich Signell wrote:

> Robb,
>
> Yes, my calcs were done using doubles (Matlab uses
> doubles by default).
>
> -Rich
>
> --- Robb Kambic <address@hidden> wrote:
> > On Thu, 1 Apr 2004, Richard Signell wrote:
> >
> > > Robb,
> > >
> > > I can't believe I didn't notice this before, but
> > you are using
> > > pole_lat, pole_lon = (-32.5, 10)
> > >
> > > but for the rtll formula I gave you to produce the
> > right results
> > > I was using
> > >
> > > pole_lat, pole_lon = (57.5, 10)
> > > which yields true lon/lat domain LL and UR
> > corners:
> > >
> > > lon(1,1)=4.5337
> > > lat(1,1)=33.2988
> > >
> > > lon(272,234)=24.6762
> > > lat(272,234)=49.4020
> >
> > Rich,
> >
> > The lats/lons are almost correct now, but I have
> > some round-off error.
> >
> > The first one is correct
> > la1 =-24, lo1 = -5
> > lon_true[ 1, 1 ] =4.53367311069887, lat_true[ 1, 1 ]
> > =33.2987847224989
> >
> > if I plug in the following values, then I get the
> > same answers you did.
> > la2 =-7.063, lo2 = 9.563
> > lon_true[ 2, 2 ] =24.67621972076, lat_true[ 2, 2 ]
> > =49.401952455077
> >
> > But iterations has round-off
> > lon[ 234 ] =10.269655932, lat[ 272 ] =-6.240013916
> > lon_true[ 234, 272 ] =26.0203076954498, lat_true[
> > 234, 272 ]
> > =50.0461679747345
> >
> > Next, I going to try doubles instead of floats.
> >
> > Robb...
> >
> >
> > >
> > > See attached image.
> > >
> > > So I guess either you need to insert a:
> > >
> > > pole_lat=90-pole_lat;
> > >
> > > line before you use my rtll formula, or else DWD
> > have
> > > incorrectly specified the "pole_lat" as
> > 90-"pole_lat" in the grib file.
> > >
> > > -Rich
> > >
> > >
> > > Robb Kambic wrote:
> > > >
> > > > On Tue, 30 Mar 2004, Richard Signell wrote:
> > > >
> > > > > Robb Kambic wrote:
> > > > > >
> > > > > > Rich,
> > > > > >
> > > > > > I have the formulas needed but the results
> > don't seem correct to me. The
> > > > > > pole is lat/lon  (-32.5, 10).  So ( -24, -5
> > ) -> ( -32.5, 10 ) and
> > > > > > (-7,063, 9.563 ) -> ( -15.563, 24.563 )
> > Since the rotation angle is 0,
> > > > > > there is no rotation.    The lat/lon
> > increment is .065535004
> > > > > >
> > > > > > The problem is the calculated corner from
> > the starting point of
> > > > > > ( -32.5, 10 ) ends up at ( -17.23, 27.76 )
> > instead of the given one
> > > > > > ( -15.563, 24.563 ).  It's a couple of
> > degrees off.  I don't think it's a
> > > > > > rounding err, so I must be missing
> > something.
> > > > > >
> > > > >
> > > > >
> > > > > I think I see your problem.   In the rotated
> > pole coordinate, the grid
> > > > > is uniformly spaced in lat/lon, but in the
> > rotated pole coordinate,
> > > > > it *is not*.   You therefore need 2D arrays to
> > hold the true lon/lat
> > > > > values, and you need to compute each of these
> > values using the formula.
> > > >
> > > > Rich,
> > > >
> > > > OK, I duplicated the code and the results still
> > are not correct. ie
> > > >
> > > > should be lon= 24.563 and lat =  -15.563
> > > >
> > > > straight increment, had error in first
> > calculation
> > > > lon[ 234 ] = 25.269655932, lat[ 272 ] =
> > -14.7400139159998
> > > >
> > > > results from duplicated rtll
> > > > almd=44.4907691654722, aphd=-43.1948333734209
> > > >
> > > > How about you send the lat lon values thru
> > MatLab and see what values you
> > > > get for the ending coordinates.
> > > >
> > > > Robb..
> > > >
> > > > >
> > > > > In other words, you need to loop thru each
> > grid point(e.g.):
> > > > >
> > > > > for j=1:272
> > > > >   lat=lat(j);
> > > > >   for i=1:234
> > > > >     lon=lon(i)
> > > > >
> >
> [lon_true(i,j),lat_true(i,j)]=rtll(pole_lon,pole_lat,lon(i),lat(i));
> > > > >   end
> > > > > end
> > > > >
> > > > > and lon_true, lat_true will be 2D arrays:
> > > > >
> > > > > lon_true(lat,lon)
> > > > > lat_true(lat,lon)
> > > > >
> > > > >
> > > > > > On Mon, 29 Mar 2004, Richard Signell wrote:
> > > > > >
> > > > > > > Robb,
> > > > > > >
> > > > > > > If it's helpful, here is my MATLAB
> > function for deriving true lon,lat
> > > > > > > from the rotated pole lon,lat.
> > > > > > >
> > > > > > > function
> > [almd,aphd]=rtll(tlm0d,tph0d,tlmd,tphd);
> > > > > > >
> >
> %-------------------------------------------------------------------------
> > > > > > > % RTLL transforms rotated coordinates
> > (tlmd,tphd) into
> > > > > > > % ordinary geographic coordinates
> > (almd,aphd). i/o decimal degrees (DD)
> > > > > > > % tlm0d, tph0d: lon e lat of the center of
> > rotation ("North Pole") in
> > > > > > > DD.
> > > > > > > %
> > > > > > > dtr=pi/180.;
> > > > > > >
> > > > > > > ctph0=cos(tph0d*dtr);
> > > > > > > stph0=sin(tph0d*dtr);
> > > > > > >
> > > > > > > stph=sin(tphd*dtr);
> > > > > > > ctph=cos(tphd*dtr);
> > > > > > > ctlm=cos(tlmd*dtr);
> > > > > > > stlm=sin(tlmd*dtr);
> > > > > > >
> > > > > > > aph=asin(stph0.*ctph.*ctlm+ctph0.*stph);
> > > > > > > cph=cos(aph);
> > > > > > >
> > > > > > > almd=tlm0d+asin(stlm.*ctph./cph)/dtr;
> > > > > > > aphd=aph/dtr;
> > > > > > > % end
> > > > > > >
> > > > > > >
> > > > > > > -Rich
> > > > > > >
> > > > > > >
> > > > >
> > > > > --
> > > > > Dr. Richard P. Signell                 |
> > address@hidden
> > > > > NATO/SACLANT Undersea Research Centre  |
> > Tel: (+39) 0187 527 381
> > > > > Viale San Bartolomeo 400               |
> > Fax: (+39) 0187 527 331
> > > > > 19138 La Spezia, ITALY  --> From USA/CANADA,
> > use: APO AE 09613-5000
> > > > >
> > > > >
> > > > > ***PRIVILEGED AND CONFIDENTIAL***
> > > > > The information contained in this e-mail
> > message (including any attached files) is intended
> > for the use of the addressee(s) only and is
> > priviliged information. The information should
> > neither be posted to the Internet, nor published in
> > any other public domain, without the express
> > permission
> === message truncated ===
>
>
> =====
> Rich Signell
> work: address@hidden
> home: address@hidden
>

===============================================================================
Robb Kambic                                Unidata Program Center
Software Engineer III                      Univ. Corp for Atmospheric Research
address@hidden             WWW: http://www.unidata.ucar.edu/
===============================================================================