COMPGROUPS.NET | Search | Post Question | Groups | Stream | About | Register

### TriScatteredInterp and duplicate data

• Email
• Follow

```Hi,

Does anyone how matlab does the average in the duplicate input of TriScatteredInterp ?
I have input a time series of scattered 2D data, so, each 2D data of this time series has its coordinates almost identical from one time to another, and a bit different values for these same coordinates. I have a warning message saying it the duplicate data were averaged. What would be the average ? Arithmetic average ? or else ?

Thanks
```
 0

See related articles to this posting

```Hi Raphael,

Yes, the average is the arithmetic average; you can verify this with a
simple example;

x = rand(100,1)*4-2;
y = rand(100,1)*4-2;
v = x.^2 + y.^2;

% Deliberately introduce duplicates so that points 25 and 50 are duplicates
of the first point.
x(25) = x(1);
y(25) = y(1);

x(50) = x(1);
y(50) = y(1);

F = TriScatteredInterp(x,y,v)

(v(1)+v(25)+v(50))/3
% The first point is the "master" point, 25 and 50 are duplicates
F.V(1)

Is this what you expected TriScatteredInterp to do, if not what did you
expect?

If you wish to develop your own custom averaging scheme you can call UNIQUE
to identify the duplicates and preprocess them before calling
TriScatteredInterp.

Best regards,

Damian

"Raphael Attie" <attie@mps.mpg.de> wrote in message
news:he3npf\$ncl\$1@fred.mathworks.com...
> Hi,
>
> Does anyone how matlab does the average in the duplicate input of
> TriScatteredInterp ?
> I have input a time series of scattered 2D data, so, each 2D data of this
> time series has its coordinates almost identical from one time to another,
> and a bit different values for these same coordinates. I have a warning
> message saying it the duplicate data were averaged. What would be the
> average ? Arithmetic average ? or else ?
>
> Thanks

```
 0

```Hi Damian,

Thanks for the prompt reply. Yes, it is indeed what I expected. It's ok like this, I just wanted to make sure I understood the warning issued by this procedure. For the duplicate, I prefer averaging myself if needed, for more transparency with my data.

Raphael

"Damian Sheehy" <Damian.Sheehy@mathworks.com> wrote in message <he3pt2\$77u\$1@fred.mathworks.com>...
> Hi Raphael,
>
> Yes, the average is the arithmetic average; you can verify this with a
> simple example;
>
> x = rand(100,1)*4-2;
> y = rand(100,1)*4-2;
> v = x.^2 + y.^2;
>
> % Deliberately introduce duplicates so that points 25 and 50 are duplicates
> of the first point.
> x(25) = x(1);
> y(25) = y(1);
>
> x(50) = x(1);
> y(50) = y(1);
>
> F = TriScatteredInterp(x,y,v)
>
> (v(1)+v(25)+v(50))/3
> % The first point is the "master" point, 25 and 50 are duplicates
> F.V(1)
>
> Is this what you expected TriScatteredInterp to do, if not what did you
> expect?
>
> If you wish to develop your own custom averaging scheme you can call UNIQUE
> to identify the duplicates and preprocess them before calling
> TriScatteredInterp.
>
> Best regards,
>
> Damian
>
>
> "Raphael Attie" <attie@mps.mpg.de> wrote in message
> news:he3npf\$ncl\$1@fred.mathworks.com...
> > Hi,
> >
> > Does anyone how matlab does the average in the duplicate input of
> > TriScatteredInterp ?
> > I have input a time series of scattered 2D data, so, each 2D data of this
> > time series has its coordinates almost identical from one time to another,
> > and a bit different values for these same coordinates. I have a warning
> > message saying it the duplicate data were averaged. What would be the
> > average ? Arithmetic average ? or else ?
> >
> > Thanks
>
```
 0

```Ooops, sent it too quickly, forgot somethin'

Damian, wouldn't be actually the use of "intersect" better to identify duplicates ? Unique would make me get the complementary set.

> "Damian Sheehy" <Damian.Sheehy@mathworks.com> wrote in message <he3pt2\$77u\$1@fred.mathworks.com>...
> > Hi Raphael,
> >
> > Yes, the average is the arithmetic average; you can verify this with a
> > simple example;
> >
> > x = rand(100,1)*4-2;
> > y = rand(100,1)*4-2;
> > v = x.^2 + y.^2;
> >
> > % Deliberately introduce duplicates so that points 25 and 50 are duplicates
> > of the first point.
> > x(25) = x(1);
> > y(25) = y(1);
> >
> > x(50) = x(1);
> > y(50) = y(1);
> >
> > F = TriScatteredInterp(x,y,v)
> >
> > (v(1)+v(25)+v(50))/3
> > % The first point is the "master" point, 25 and 50 are duplicates
> > F.V(1)
> >
> > Is this what you expected TriScatteredInterp to do, if not what did you
> > expect?
> >
> > If you wish to develop your own custom averaging scheme you can call UNIQUE
> > to identify the duplicates and preprocess them before calling
> > TriScatteredInterp.
> >
> > Best regards,
> >
> > Damian
> >
> >
> > "Raphael Attie" <attie@mps.mpg.de> wrote in message
> > news:he3npf\$ncl\$1@fred.mathworks.com...
> > > Hi,
> > >
> > > Does anyone how matlab does the average in the duplicate input of
> > > TriScatteredInterp ?
> > > I have input a time series of scattered 2D data, so, each 2D data of this
> > > time series has its coordinates almost identical from one time to another,
> > > and a bit different values for these same coordinates. I have a warning
> > > message saying it the duplicate data were averaged. What would be the
> > > average ? Arithmetic average ? or else ?
> > >
> > > Thanks
> >
```
 0

```Hi Raphael,

I am not sure I fully understand your usecase. There are two scenarios where
you are likely to encounter duplicate points; one is during the construction
of the interpolant and I just provided an example of that case, the other
case may arise when you insert additional datapoints into an existing
interpolant.

I am assuming you are dealing with the first scenario and in this case I
believe the function UNIQUE would be a suitable choice.

If you are dealing with the second scenario then the INTERSECT function
would be the one to use. If this is indeed your usecase I would be inclined
to work directly with DelaunayTri IF you happen to be doing "nearest" or
"linear" interpolation. It would be a few extra steps, but you may be able
to get a more efficient solution if performance is important. Let me know if
you would like the details.

Regards,

Damian

"Raphael Attie" <attie@mps.mpg.de> wrote in message
news:he3sb7\$ar2\$1@fred.mathworks.com...
> Ooops, sent it too quickly, forgot somethin'
>
> Damian, wouldn't be actually the use of "intersect" better to identify
> duplicates ? Unique would make me get the complementary set.
>
>
>> "Damian Sheehy" <Damian.Sheehy@mathworks.com> wrote in message
>> <he3pt2\$77u\$1@fred.mathworks.com>...
>> > Hi Raphael,
>> >
>> > Yes, the average is the arithmetic average; you can verify this with a
>> > simple example;
>> >
>> > x = rand(100,1)*4-2;
>> > y = rand(100,1)*4-2;
>> > v = x.^2 + y.^2;
>> >
>> > % Deliberately introduce duplicates so that points 25 and 50 are
>> > duplicates
>> > of the first point.
>> > x(25) = x(1);
>> > y(25) = y(1);
>> >
>> > x(50) = x(1);
>> > y(50) = y(1);
>> >
>> > F = TriScatteredInterp(x,y,v)
>> >
>> > (v(1)+v(25)+v(50))/3
>> > % The first point is the "master" point, 25 and 50 are duplicates
>> > F.V(1)
>> >
>> > Is this what you expected TriScatteredInterp to do, if not what did you
>> > expect?
>> >
>> > If you wish to develop your own custom averaging scheme you can call
>> > UNIQUE
>> > to identify the duplicates and preprocess them before calling
>> > TriScatteredInterp.
>> >
>> > Best regards,
>> >
>> > Damian
>> >
>> >
>> > "Raphael Attie" <attie@mps.mpg.de> wrote in message
>> > news:he3npf\$ncl\$1@fred.mathworks.com...
>> > > Hi,
>> > >
>> > > Does anyone how matlab does the average in the duplicate input of
>> > > TriScatteredInterp ?
>> > > I have input a time series of scattered 2D data, so, each 2D data of
>> > > this
>> > > time series has its coordinates almost identical from one time to
>> > > another,
>> > > and a bit different values for these same coordinates. I have a
>> > > warning
>> > > message saying it the duplicate data were averaged. What would be the
>> > > average ? Arithmetic average ? or else ?
>> > >
>> > > Thanks
>> >

```
 0

```Hi there,

Indeed, performance is an issue.

As for the unique/intersect problem: well, i have velocity fields. Each field given at a specific time, on its own irregular grid (picked up by some semi-Lagrangien's probes that move in time).  This irregular grid varies a little in time, so, lots of coordinates are identical from one time to another and my goal is to have one final averaged uniformly gridded velocity field, that would be the average of this time series. I've already done it with the routine, but it's quite slow, my velocity fields are meant to be 1004x498 each, and I have about 100 of them. I was hoping about something faster, by avoiding this average done within the TriScatteredInterp...

Do you think that doing things separately,  with the separate Delaunay triangulation, would make things faster ? What do you think would be the best approach ?
Other alternative : is TriScatteredInterp using the multi-cpus ? If I run it on bigger machines at my working place, I'd be using up to 12 cpus, would it speed things up ? Right now i'm testing on dual-core machine only.

Thanks again, it's very nice to have someonelse's opinion.

Raphael

"Damian Sheehy" <Damian.Sheehy@mathworks.com> wrote in message <he4aa4\$erq\$1@fred.mathworks.com>...
> Hi Raphael,
>
> I am not sure I fully understand your usecase. There are two scenarios where
> you are likely to encounter duplicate points; one is during the construction
> of the interpolant and I just provided an example of that case, the other
> case may arise when you insert additional datapoints into an existing
> interpolant.
>
> I am assuming you are dealing with the first scenario and in this case I
> believe the function UNIQUE would be a suitable choice.
>
> If you are dealing with the second scenario then the INTERSECT function
> would be the one to use. If this is indeed your usecase I would be inclined
> to work directly with DelaunayTri IF you happen to be doing "nearest" or
> "linear" interpolation. It would be a few extra steps, but you may be able
> to get a more efficient solution if performance is important. Let me know if
> you would like the details.
>
>
> Regards,
>
> Damian
>
>
> "Raphael Attie" <attie@mps.mpg.de> wrote in message
> news:he3sb7\$ar2\$1@fred.mathworks.com...
> > Ooops, sent it too quickly, forgot somethin'
> >
> > Damian, wouldn't be actually the use of "intersect" better to identify
> > duplicates ? Unique would make me get the complementary set.
> >
> >
> >> "Damian Sheehy" <Damian.Sheehy@mathworks.com> wrote in message
> >> <he3pt2\$77u\$1@fred.mathworks.com>...
> >> > Hi Raphael,
> >> >
> >> > Yes, the average is the arithmetic average; you can verify this with a
> >> > simple example;
> >> >
> >> > x = rand(100,1)*4-2;
> >> > y = rand(100,1)*4-2;
> >> > v = x.^2 + y.^2;
> >> >
> >> > % Deliberately introduce duplicates so that points 25 and 50 are
> >> > duplicates
> >> > of the first point.
> >> > x(25) = x(1);
> >> > y(25) = y(1);
> >> >
> >> > x(50) = x(1);
> >> > y(50) = y(1);
> >> >
> >> > F = TriScatteredInterp(x,y,v)
> >> >
> >> > (v(1)+v(25)+v(50))/3
> >> > % The first point is the "master" point, 25 and 50 are duplicates
> >> > F.V(1)
> >> >
> >> > Is this what you expected TriScatteredInterp to do, if not what did you
> >> > expect?
> >> >
> >> > If you wish to develop your own custom averaging scheme you can call
> >> > UNIQUE
> >> > to identify the duplicates and preprocess them before calling
> >> > TriScatteredInterp.
> >> >
> >> > Best regards,
> >> >
> >> > Damian
> >> >
> >> >
> >> > "Raphael Attie" <attie@mps.mpg.de> wrote in message
> >> > news:he3npf\$ncl\$1@fred.mathworks.com...
> >> > > Hi,
> >> > >
> >> > > Does anyone how matlab does the average in the duplicate input of
> >> > > TriScatteredInterp ?
> >> > > I have input a time series of scattered 2D data, so, each 2D data of
> >> > > this
> >> > > time series has its coordinates almost identical from one time to
> >> > > another,
> >> > > and a bit different values for these same coordinates. I have a
> >> > > warning
> >> > > message saying it the duplicate data were averaged. What would be the
> >> > > average ? Arithmetic average ? or else ?
> >> > >
> >> > > Thanks
> >> >
>
```
 0

```Ohhh, OK, I get it, you have something like an ALE moving mesh simulation.
Yes, I think you will be far better off working directly with DelaunayTri,
but you will be limited to linear interpolation. Unfortunately you will not
be able to leverage multicore, but you can still do better.

TriScatteredInterp uses an underlying Delaunay triangulation to perform the
interpolation.
The bottleneck in your approach is that you are recreating the entire
triangulation during each timestep, wereas only a fraction of the nodes have
moved and the triangulation only needs to be updated locally.

Here's what you can do. . .

First understand how to perform interpolation using the DelaunayTri.
Lets say your sample data is x, y, v as in the case of TriScatteredInterp
and the grid points where you want to interpolate are given by xi, yi.

% Create the Delaunay
dt = DelaunayTri(x,y)

% To Interpolate at a *single* query points xiq, yq you do the following

[T, BC] = dt.pointLocation(xq,yq)

% This locates the indices of the triangle containing (xq, yq) and the
% corresponding Barycentric Coordinates (weights associated with each point)

% The value at xq, yq is given by the
% Value-at-first-vertex*the-first-weight +
v(dt(T,1))*BC(1) + v(dt(T,2))*BC(2) + v(dt(T,3))*BC(3)

% This can be vectorized for a set of query points (xi,yi)
% You should check this against what you get from TriScatteredInterp
[T, BC] = dt.pointLocation(xi,yi)
vi = v(dt(T,1)).*BC(:,1) +  v(dt(T,2)).*BC(:,2) +  v(dt(T,3)).*BC(:,3)

OK, next step is to update the triangulation to capture the nodes that moved
to xt2, yt2
Ignore the nodes by moved by such a small amount (DistTolerance) that the
movement can be ignored.

pid = dt.nearestNeighbor(xt2,yt2);
% Compute the distance from (xt2,yt2) to dt.X(pid,:) and (Yes,
nearestNeighbor should have have optionally provided this. I will fix it...)

If nodes J moved a distance >  DistTolerance update their new locations

dt.X(J,:) = [xt2(J),yt2(J)]

Now interpolate at t2, just as before
[T, BC] = dt.pointLocation(xi,yi)
vit2 = v(dt(T,1)).*BC(:,1) +  v(dt(T,2)).*BC(:,2) +  v(dt(T,3)).*BC(:,3)

If your finite elememt mesh is composed of triangular elements you could
interpolate on that also but it's a bit more involved.

Send me an email if you have any further questions. . .

All the best,

Damian

"Raphael Attie" <attie@mps.mpg.de> wrote in message
news:he4qff\$66f\$1@fred.mathworks.com...
> Hi there,
>
> Indeed, performance is an issue.
>
> As for the unique/intersect problem: well, i have velocity fields. Each
> field given at a specific time, on its own irregular grid (picked up by
> some semi-Lagrangien's probes that move in time).  This irregular grid
> varies a little in time, so, lots of coordinates are identical from one
> time to another and my goal is to have one final averaged uniformly
> gridded velocity field, that would be the average of this time series.
> I've already done it with the routine, but it's quite slow, my velocity
> fields are meant to be 1004x498 each, and I have about 100 of them. I was
> hoping about something faster, by avoiding this average done within the
> TriScatteredInterp...
>
> Do you think that doing things separately,  with the separate Delaunay
> triangulation, would make things faster ? What do you think would be the
> best approach ?
> Other alternative : is TriScatteredInterp using the multi-cpus ? If I run
> it on bigger machines at my working place, I'd be using up to 12 cpus,
> would it speed things up ? Right now i'm testing on dual-core machine
> only.
>
> Thanks again, it's very nice to have someonelse's opinion.
>
> Raphael
>
> "Damian Sheehy" <Damian.Sheehy@mathworks.com> wrote in message
> <he4aa4\$erq\$1@fred.mathworks.com>...
>> Hi Raphael,
>>
>> I am not sure I fully understand your usecase. There are two scenarios
>> where
>> you are likely to encounter duplicate points; one is during the
>> construction
>> of the interpolant and I just provided an example of that case, the other
>> case may arise when you insert additional datapoints into an existing
>> interpolant.
>>
>> I am assuming you are dealing with the first scenario and in this case I
>> believe the function UNIQUE would be a suitable choice.
>>
>> If you are dealing with the second scenario then the INTERSECT function
>> would be the one to use. If this is indeed your usecase I would be
>> inclined
>> to work directly with DelaunayTri IF you happen to be doing "nearest" or
>> "linear" interpolation. It would be a few extra steps, but you may be
>> able
>> to get a more efficient solution if performance is important. Let me know
>> if
>> you would like the details.
>>
>>
>> Regards,
>>
>> Damian
>>
>>
>> "Raphael Attie" <attie@mps.mpg.de> wrote in message
>> news:he3sb7\$ar2\$1@fred.mathworks.com...
>> > Ooops, sent it too quickly, forgot somethin'
>> >
>> > Damian, wouldn't be actually the use of "intersect" better to identify
>> > duplicates ? Unique would make me get the complementary set.
>> >
>> >
>> >> "Damian Sheehy" <Damian.Sheehy@mathworks.com> wrote in message
>> >> <he3pt2\$77u\$1@fred.mathworks.com>...
>> >> > Hi Raphael,
>> >> >
>> >> > Yes, the average is the arithmetic average; you can verify this with
>> >> > a
>> >> > simple example;
>> >> >
>> >> > x = rand(100,1)*4-2;
>> >> > y = rand(100,1)*4-2;
>> >> > v = x.^2 + y.^2;
>> >> >
>> >> > % Deliberately introduce duplicates so that points 25 and 50 are
>> >> > duplicates
>> >> > of the first point.
>> >> > x(25) = x(1);
>> >> > y(25) = y(1);
>> >> >
>> >> > x(50) = x(1);
>> >> > y(50) = y(1);
>> >> >
>> >> > F = TriScatteredInterp(x,y,v)
>> >> >
>> >> > (v(1)+v(25)+v(50))/3
>> >> > % The first point is the "master" point, 25 and 50 are duplicates
>> >> > F.V(1)
>> >> >
>> >> > Is this what you expected TriScatteredInterp to do, if not what did
>> >> > you
>> >> > expect?
>> >> >
>> >> > If you wish to develop your own custom averaging scheme you can call
>> >> > UNIQUE
>> >> > to identify the duplicates and preprocess them before calling
>> >> > TriScatteredInterp.
>> >> >
>> >> > Best regards,
>> >> >
>> >> > Damian
>> >> >
>> >> >
>> >> > "Raphael Attie" <attie@mps.mpg.de> wrote in message
>> >> > news:he3npf\$ncl\$1@fred.mathworks.com...
>> >> > > Hi,
>> >> > >
>> >> > > Does anyone how matlab does the average in the duplicate input of
>> >> > > TriScatteredInterp ?
>> >> > > I have input a time series of scattered 2D data, so, each 2D data
>> >> > > of
>> >> > > this
>> >> > > time series has its coordinates almost identical from one time to
>> >> > > another,
>> >> > > and a bit different values for these same coordinates. I have a
>> >> > > warning
>> >> > > message saying it the duplicate data were averaged. What would be
>> >> > > the
>> >> > > average ? Arithmetic average ? or else ?
>> >> > >
>> >> > > Thanks
>> >> >
>>

```
 0

```Hi again Damian,

I haven't been to my workstation yet but I wonder why the griddata function is not recommended any more ? As it is also based on CGAL. For what griddata was doing, why isn't it recommended anymore ? Is this just efficiency issues ? are the results the same ? if not, then why not ?

Thanks

Raphael
```
 0

```Hi  Damian,

I too work on 3D unstructured grid (tetrahedral) generated from COMSOL/FEM.
I saw your explanation on TriScatterdInterp using the DelaunayTri.
I really appreciate your understanding. It would also be nice if you can explain the 'natural'
choice of interpolation as you did for 'nearest' and 'linear' using the [T , BC] ...

Can you also suggest how to extend it to 'bicubic' or any higher order method. Because current  "TriScatterdInterp" does only linear,nearest and natural.

Thanks
Neehar
```
 0