[MITgcm-support] tsearch

[Yuan Lian]

Hi Martin,

Here is the version that may be checked into the repository. I also 
added a version check in case someone is using a older version of 
matlab. Please note that I assumed version 7.14 or later removed tsearch.

Best,
Yuan


On 7/9/13 12:20 AM, Martin Losch wrote:
> Hi Angela,
>
> could you send me a working version of the matlab-script (preferably with the not-working code lines commented out followed by the new code lines)? I could then check this into the repository.
>
> Martin
>
> On Jul 9, 2013, at 2:11 AM, Angela Zalucha <azalucha at seti.org> wrote:
>
>> Thanks, that did work.
>>
>>   Angela
>>
>> _-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
>> Angela Zalucha, Ph.D.
>> Research Scientist
>> Search for Extraterrestrial Intelligence (SETI) Institute
>>
>>
>> Office located at Southwest Research Institute
>> 1050 Walnut Street, Suite 300
>> Boulder, CO 80302
>> USA
>> (720) 208-7211
>>
>>
>> http://www.seti.org/users/azalucha
>> _-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
>>
>> On Mon, 8 Jul 2013, Yuan Lian wrote:
>>
>>> Hi Angela,
>>>
>>> You can replace tsearch line with the following:
>>>
>>> %%%%%% the following only works for matlab version 2012 or later.
>>> % Triangularize the data
>>> tri=DelaunayTri(x,y);
>>> if isempty(tri),
>>> warning('Data cannot be triangulated.');
>>> return
>>> end
>>>
>>> % Find the nearest triangle (t)
>>> [t, bcs] = pointLocation(tri, xi, yi);
>>>
>>> %%%%%%
>>>
>>> Hope this helps.
>>>
>>> Best,
>>> Yuan
>>>
>>>
>>> On 7/8/13 3:12 PM, Angela Zalucha wrote:
>>>> Yes, I have found this page:
>>>> http://www.mathworks.com/matlabcentral/fileexchange/authors/251793
>>>> But I don't know what "p" corresponds to in griddata_preprocess.m, as it is a black box to me.
>>>>
>>>>   Angela
>>>> _-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
>>>> Angela Zalucha, Ph.D.
>>>> Research Scientist
>>>> Search for Extraterrestrial Intelligence (SETI) Institute
>>>> Office located at Southwest Research Institute
>>>> 1050 Walnut Street, Suite 300
>>>> Boulder, CO 80302
>>>> USA
>>>> (720) 208-7211
>>>> http://www.seti.org/users/azalucha
>>>> _-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
>>>> On Mon, 8 Jul 2013, Menemenlis, Dimitris (3248) wrote:
>>>>> Angela you are right, tsearch has been replaced by DelaunayTri function in recent matlab.
>>>>> You could try replacing tsearch with DelaunayTri function in griddata_preprocess.m
>>>>> But this routine is not needed for reading MITgcm files.
>>>>> It is used for computing interpolation weights.
>>>>> Dimitris Menemenlis
>>>>> On Jul 8, 2013, at 5:50 PM, Angela Zalucha wrote:
>>>>>> Hi,
>>>>>> I am running Matlab 2013a and I am unable to read in any files due to a Matlab function that is no longer supported, tsearch.  This function is called in the mitgcm-provided script griddata_preprocess.m.  I checked the repository and there is no update there.  Has anyone fixed this problem?
>>>>>>
>>>>>> Angela
>>>>> _______________________________________________
>>>>> MITgcm-support mailing list
>>>>> MITgcm-support at mitgcm.org
>>>>> http://mitgcm.org/mailman/listinfo/mitgcm-support
>>>> _______________________________________________
>>>> MITgcm-support mailing list
>>>> MITgcm-support at mitgcm.org
>>>> http://mitgcm.org/mailman/listinfo/mitgcm-support
>>>
>>> _______________________________________________
>>> MITgcm-support mailing list
>>> MITgcm-support at mitgcm.org
>>> http://mitgcm.org/mailman/listinfo/mitgcm-support
>>>
>> _______________________________________________
>> MITgcm-support mailing list
>> MITgcm-support at mitgcm.org
>> http://mitgcm.org/mailman/listinfo/mitgcm-support
>
> _______________________________________________
> MITgcm-support mailing list
> MITgcm-support at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-support

-------------- next part --------------
function [del] = griddata_preprocess(x,y,xi,yi,method)
%GRIDDATA_PREPROCESS Pre-calculate Delaunay triangulation for use
%   with GRIDDATA_FAST.
%
%   DEL = GRIDDATA_PREPROCESS(X,Y,XI,YI)

%   Clay M. Thompson 8-21-95
%   Copyright 1984-2001 The MathWorks, Inc. 
%   $Revision: 1.1.1.1 $  $Date: 2003/11/11 18:08:08 $

error(nargchk(4,5,nargin))

if prod(size(xi)) ~= prod(size(yi))
 [yi,xi]=ndgrid(yi,xi);
end

if nargin<6, method = 'linear'; end
if ~isstr(method), 
  error('METHOD must be one of ''linear'',''cubic'',''nearest'', or ''v4''.');
end


switch lower(method),
  case 'linear'
    del = linear(x,y,xi,yi);
% case 'cubic'
%   zi = cubic(x,y,z,xi,yi);
% case 'nearest'
%   zi = nearest(x,y,z,xi,yi);
% case {'invdist','v4'}
%   zi = gdatav4(x,y,z,xi,yi);
  otherwise
    error('Unknown method.');
end
  


%------------------------------------------------------------
function delau = linear(x,y,xi,yi)
%LINEAR Triangle-based linear interpolation

%   Reference: David F. Watson, "Contouring: A guide
%   to the analysis and display of spacial data", Pergamon, 1994.

siz = size(xi);
xi = xi(:); yi = yi(:); % Treat these as columns
x = x(:); y = y(:); % Treat these as columns

%%%%%% the following only works for matlab version 2011 or older.
if verLessThan('matlab','7.14') % 7.14 is 2012a ???

   % Triangularize the data
   tri = delaunayn([x y]);
   if isempty(tri),
      warning('Data cannot be triangulated.');
      return
   end

   % Find the nearest triangle (t)
   t = tsearch(x,y,tri,xi,yi);

%%%%%% the following only works for matlab version 2012 or later.
else % version is 2012 or newer

   % Triangularize the data
   tri=DelaunayTri(x,y);
   if isempty(tri),
      warning('Data cannot be triangulated.');
      return
   end

   % Find the nearest triangle (t)
   [t, bcs] = pointLocation(tri, xi, yi);

%%%%%%
end; % end of selecting version

% Only keep the relevant triangles.
out = find(isnan(t));
if ~isempty(out), t(out) = ones(size(out)); end
tri = tri(t,:);

% Compute Barycentric coordinates (w).  P. 78 in Watson.
del = (x(tri(:,2))-x(tri(:,1))) .* (y(tri(:,3))-y(tri(:,1))) - ...
      (x(tri(:,3))-x(tri(:,1))) .* (y(tri(:,2))-y(tri(:,1)));
w(:,3) = ((x(tri(:,1))-xi).*(y(tri(:,2))-yi) - ...
          (x(tri(:,2))-xi).*(y(tri(:,1))-yi)) ./ del;
w(:,2) = ((x(tri(:,3))-xi).*(y(tri(:,1))-yi) - ...
          (x(tri(:,1))-xi).*(y(tri(:,3))-yi)) ./ del;
w(:,1) = ((x(tri(:,2))-xi).*(y(tri(:,3))-yi) - ...
          (x(tri(:,3))-xi).*(y(tri(:,2))-yi)) ./ del;
w(out,:) = zeros(length(out),3);

delau.tri=tri;
delau.w=w;
delau.siz=siz;
delau.out=out;

%------------------------------------------------------------