Search All of the Math Forum:
Views expressed in these public forums are not endorsed by
NCTM or The Math Forum.


rachael
Posts:
16
Registered:
2/10/09


Re: griddata error: qhull precision error, ERRONEOUS FACET
Posted:
Feb 24, 2009 3:52 PM


"rachael " <rmueller@coas.oregonstate.edu> wrote in message <gnso7i$k98$1@fred.mathworks.com>... > I am trying to grid data 1 arc second grid onto a 4 Km grid, using griddata. I have tried inputting the data twice, and I get a similar error both times. The first time, I used a version of the data the comes from someone else; the second time, I am using a version of the data that I created from a downloaded dataset. The data that I "created" is from the GEBCO dataset of gridded bathymetry. > <http://www.gebco.net/data_and_products/gridded_bathymetry_data/> > It's 1 arcminute, with equidistant grid spacing. I create my latitude (LAT), longitude(LON) and bathymetry (Z) grids as follows: > > > ncload('/Users/rmueller/Data/GEBCO /62s75s80w40w_GEBCO_1min.nc','x_range','y_range','z_range','z'); > > > ncload('/Users/rmueller/Data/GEBCO/62s75s80w40w_GEBCO_1min.nc','spacing','dimension'); > > > dx = spacing(1);dy=spacing(2); > > > x=[x_range(2):1*dx:x_range(1)];y=[y_range(2):1*dy:y_range(1)]; > > > [LAT,LON]=meshgrid(y,x); > > > Z = flipud(reshape(z,dimension(1),dimension(2))); > > The errors arise when I use griddata to map this data onto a 4km spaced, polar stereographic, lat_rho/lon_rho grid using the following > >z = griddata(LON,LAT,Z,lon_rho,lat_rho); > > I have used this lat_rho, lon_rho to convert a different dataset w/o problems, so I think it's the GEBCO data. Since I get a similar error with the data as converted into matlab format by two different users, I think the problem is more due to the data than the format I have it in matlab. But what is the problem??? Is the 1 arcminute too great of a grid spacing to map onto a 4 Km grid??? > > The two errors are shown in detail below. > > ERROR 1: > ??? Error using ==> qhullmx > qhull precision error: initial facet 2 is coplanar with the interior point > > ERRONEOUS FACET: > While executing:  qhull d Qt Qbb Qc > Options selected for Qhull 2003.1 2003/12/30: > delaunay Qtriangulate Qbboundlast Qcoplanarkeep _premerge > _zerocentrum Pgood Qinteriorkeep _maxwidth 35 Errorroundoff 1.1e13 > _onemerge 7.4e13 Visibledistance 2.1e13 Ucoplanardistance 2.1e13 > Widthoutside 4.2e13 _widefacet 1.3e12 > precision problems (corrected unless 'Q0' or an error) > 2 > flipped facets > The input to qhull appears to be less than 3 dimensional, or a > computation has overflowed. > Qhull could not construct a clearly convex simplex from points: > The center point is coplanar with a facet, or a vertex is coplanar > with a neighboring facet. The maximum round off error for > computing distances is 1.1e13. The center point, facets and distances > to the center point are as follows: > > facet > p2016960 > p2100 > p0 > distance= 3.6e15 > > facet > p2019060 > p2100 > p0 > distance= 0 > > facet > p2019060 > p2016960 > p0 > distance= 1.4e14 > > facet > p2019060 > p2016960 > p2100 > distance= 5.3e15 > > These points either have a maximum or minimum xcoordinate, or > they maximize the determinant for k coordinates. Trial points > are first selected from points that maximize a coordinate. > The min and max coordinates for each dimension are: > 0: 75 2.225e308 difference= 75 > 1: 76 2.225e308 difference= 76 > 2: 3.553e15 35 difference= 35 > If the input should be full dimensional, you have several options that > may determine an initial simplex: >  use 'QJ' to joggle the input and make it full dimensional >  use 'QbB' to scale the points to the unit cube >  use 'QR0' to randomly rotate the input for different maximum points >  use 'Qs' to search all points for the initial simplex >  use 'En' to specify a maximum roundoff error less than 1.1e13. >  trace execution with 'T3' to see the determinant for each point. > If the input is lower dimensional: >  use 'QJ' to joggle the input and make it full dimensional >  use 'Qbk:0Bk:0' to delete coordinate k from the input. You should > pick the coordinate with the least range. The hull will have the > correct topology. >  determine the flat containing the points, rotate the points > into a coordinate plane, and delete the other coordinates. >  add one or more points to make the input full dimensional. > This is a Delaunay triangulation and the input is cocircular or cospherical: >  use 'Qz' to add a point "at infinity" (i.e., above the paraboloid) >  or use 'QJ' to joggle the input and avoid cocircular data > Error in ==> delaunayn at 117 > t = qhullmx(x', 'd ', opt); > Error in ==> griddata>linear at 151 > tri = delaunayn([x y]); > Error in ==> griddata at 120 > zi = linear(x,y,z,xi,yi,opt); > > > ERROR 2: > ??? Error using ==> qhullmx > While executing:  qhull d Qt Qbb Qc > Options selected for Qhull 2003.1 2003/12/30: > delaunay Qtriangulate Qbboundlast Qcoplanarkeep _premerge > _zerocentrum Pgood Qinteriorkeep _maxwidth 40 Errorroundoff 1.1e13 > _onemerge 7.8e13 Visibledistance 2.2e13 Ucoplanardistance 2.2e13 > Widthoutside 4.4e13 _widefacet 1.3e12 > Last point added to hull was p383250. > Last merge was #3109857. > Qhull has finished constructing the hull. > At error exit: > precision problems (corrected unless 'Q0' or an error) > 3109857 > coplanar horizon facets for new vertices > Error in ==> delaunayn at 117 > t = qhullmx(x', 'd ', opt); > Error in ==> griddata>linear at 151 > tri = delaunayn([x y]); > Error in ==> griddata at 120 > zi = linear(x,y,z,xi,yi,opt); > > Thanks.
Problem solved! Now that I know the answer, the cryptic error message makes more sense. The 'QJ' mentioned above is used as an "option" in griddata. The solution that worked for me is:
z = griddata(lon_g,lat_g,z_g,lon_rho,lat_rho, 'linear',{'QJ'});



