[SciPy-User] Creating meshgrid from meshgrid

classic Classic list List threaded Threaded
5 messages Options
Reply | Threaded
Open this post in threaded view
|

[SciPy-User] Creating meshgrid from meshgrid

Florian Lindner
Hello,

yes, the subject sound a bit weird...

I have two arrays of size N (let's say 2 arrays of length 4) that I combine using np.meshgrid

xxA, yyA = np.meshgrid(xA, yA)
xxB, yyB = np.meshgrid(xB, yB)

which gives me two meshes

xx.shape = yy.shape = (4,4)
which represent a N-dimensional mesh with 16 elements.

Now I want to evaluate a function f on every possible pair of N-dimensional points in the grid, resulting in a 16 x 16
matrix:

in a flattened notation, pA = (xxA, yyA)

f(pA[1]-pB[1]) f(pA[1]-pB[2]) f(pA[1]-pB[3]) ...
f(pA[2]-pB[1]) f(pA[2]-pB[2]) f(pA[2]-pB[3]) ...
f(pA[3]-pB[1]) f(pA[3]-pB[2]) f(pA[3]-pB[3]) ...
.
.
.


Let's say xA = yA = [1,2,3] and xB = yB = [10,20,30]

that gives me a mesh A:

(1,3) (2,3) (3,3)
(1,2) (2,2) (3,2)
(1,1) (2,1) (3,1)

and a mesh B alike.

My result matrix now should be of size 9 x 9:

f( (1,3), (10,30) ) f( (2,3), (20,30) ) f( (3,3), (30, 30) )
f( (1,2), (10,20) ) f( (2,2), (20,20) ) f( (3,2), (30, 20) )
...


f always takes two N-dimensional vectors and returns a scalar.

I hope I was able to explain what I want to achieve.

What is the best way to do that in numpy/scipy?

As long as the meshes itself are 1-d I did it like that:

mgrid = np.meshgrid([1,2,3], [10,20,30])
A = f( np.abs(mgrid[0] - mgrid[1]) )

Thanks,
Florian
_______________________________________________
SciPy-User mailing list
[hidden email]
https://mail.python.org/mailman/listinfo/scipy-user
Reply | Threaded
Open this post in threaded view
|

Re: Creating meshgrid from meshgrid

Chris Barker - NOAA Federal
I'm confused as to what you are trying to do.

Meshgrid is used to create a regular grid when you have a vectors of the axes.

It support 2 or more dimensions.

I have two arrays of size N (let's say 2 arrays of length 4) that I combine using np.meshgrid

xxA, yyA = np.meshgrid(xA, yA)
xxB, yyB = np.meshgrid(xB, yB)

which gives me two meshes

xx.shape = yy.shape = (4,4)
which represent a N-dimensional mesh with 16 elements.

no -- it represents a 2-dimensional mesh with four nodes in each direction.
 
Now I want to evaluate a function f on every possible pair of N-dimensional points in the grid, resulting in a 16 x 16
matrix:

I think you are looking for a different function than meshgrid.

But if you want to evaluate a function on a four dimensional space, you can use meshgrid with four dimensions:

xx, yy, zz, tt = meshgrid(x, y, z, t)

results = func(xx,yy,zz,tt)

note that with numpy's broadcasting, you may not need to use meshgrid at all.

Is that what you are looking for?

-CHB


--

Christopher Barker, Ph.D.
Oceanographer

Emergency Response Division
NOAA/NOS/OR&R            (206) 526-6959   voice
7600 Sand Point Way NE   (206) 526-6329   fax
Seattle, WA  98115       (206) 526-6317   main reception

[hidden email]

_______________________________________________
SciPy-User mailing list
[hidden email]
https://mail.python.org/mailman/listinfo/scipy-user
Reply | Threaded
Open this post in threaded view
|

Re: Creating meshgrid from meshgrid

Florian Lindner
Hey,

I'm sorry, I realized my posting was confusing, tried to rewrite it:


In 1-d I can easily and elegantly achieve what I want, but I'm unable to transfer it to 2-d or even higher dimensional

**1-d**

I have two 1-d meshes of size N = 3.

    a = [1, 2, 3]
    b = [10, 20, 30]

now I want to evaluate a function fun on on the difference/norm of all pairs of these meshes:

    def fun(x):
        return x

    mgrid = np.meshgrid( a,b )
    A = fun( np.abs(mgrid[0] - mgrid[1] ) )

Result A is of size N x N:

    array([[ 9,  8,  7],
           [19, 18, 17],
           [29, 28, 27]])

which is

    |b[0]-a[0]| |b[0]-a[1]| |b[0]-a[2]|
    |b[1]-a[0]| |b[1]-a[1]| |b[1]-a[2]|
    |b[2]-a[0]| |b[2]-a[1]| |b[2]-a[2]|

(only arguments of fun, for sake of brevity)

**2-d**

Now in 2-d function fun stays the same.

I have again two meshes a and b

    ax = [1, 2, 3]
    ay = [4, 5, 6]
    bx = [10, 20, 30]
    by = [40, 50, 60]

    a = np.meshgrid(ax, ay)
    b = np.meshgrid(bx, by)

Now, how can I do the same I did with the 1-d meshes above for 2-d and possibly also for higher dimensions?

The first row of A:

    || (10 40) - (1 4) || || (10 40) - (1 5) || || (10 40) - (1 6) ||
    || (10 40) - (2 4) || || (10 40) - (2 5) || || (10 40) - (2 6) ||
    || (10 40) - (3 4) || || (10 40) - (3 5) || || (10 40) - (3 6) ||

(everything is jus the first row of A, again only arguments of fun)

The result mesh should have the size N * N x N * N.

I tried to create the coordinates of a and b

    ca = np.array(list(zip(a[0].flatten(), a[1].flatten())))
    cb = np.array(list(zip(b[0].flatten(), b[1].flatten())))

And create a meshgrid from that:

    mgrid = np.meshgrid([ca], [cb])

but alone the dimensionality does not fit (18 instead of 9).

I hope I was now able to better get across what I want.

Thanks!
Florian


Am 06.09.2017 um 07:54 schrieb Chris Barker:

> I'm confused as to what you are trying to do.
>
> Meshgrid is used to create a regular grid when you have a vectors of the axes.
>
> It support 2 or more dimensions.
>
>     I have two arrays of size N (let's say 2 arrays of length 4) that I combine using np.meshgrid
>
>     xxA, yyA = np.meshgrid(xA, yA)
>     xxB, yyB = np.meshgrid(xB, yB)
>
>     which gives me two meshes
>
>     xx.shape = yy.shape = (4,4)
>     which represent a N-dimensional mesh with 16 elements.
>
>
> no -- it represents a 2-dimensional mesh with four nodes in each direction.
>  
>
>     Now I want to evaluate a function f on every possible pair of N-dimensional points in the grid, resulting in a 16 x 16
>     matrix:
>
>
> I think you are looking for a different function than meshgrid.
>
> But if you want to evaluate a function on a four dimensional space, you can use meshgrid with four dimensions:
>
> xx, yy, zz, tt = meshgrid(x, y, z, t)
>
> results = func(xx,yy,zz,tt)
>
> note that with numpy's broadcasting, you may not need to use meshgrid at all.
>
> Is that what you are looking for?
>
> -CHB
>
>
> --
>
> Christopher Barker, Ph.D.
> Oceanographer
>
> Emergency Response Division
> NOAA/NOS/OR&R            (206) 526-6959   voice
> 7600 Sand Point Way NE   (206) 526-6329   fax
> Seattle, WA  98115       (206) 526-6317   main reception
>
> [hidden email] <mailto:[hidden email]>
>
>
> _______________________________________________
> SciPy-User mailing list
> [hidden email]
> https://mail.python.org/mailman/listinfo/scipy-user
>
_______________________________________________
SciPy-User mailing list
[hidden email]
https://mail.python.org/mailman/listinfo/scipy-user
Reply | Threaded
Open this post in threaded view
|

Re: Creating meshgrid from meshgrid

Pauli Virtanen-3
Florian Lindner kirjoitti 06.09.2017 klo 06:31:
[clip]

>      def fun(x):
>          return x
>
>      mgrid = np.meshgrid( a,b )
>      A = fun( np.abs(mgrid[0] - mgrid[1] ) )
>
> Result A is of size N x N:
>
>      array([[ 9,  8,  7],
>             [19, 18, 17],
>             [29, 28, 27]])

Don't use meshgrid, use broadcasting:

from numpy import newaxis

ax = np.array([0, 1, 2])
ay = np.array([4, 5])
bx = np.array([7, 8, 9, 10])
by = np.array([11, 12, 13, 14])

def fun(x1, y1, x2, y2):
     return abs(x1 - x2) + abs(y1 - y2)

d = fun(ax[:,newaxis,newaxis,newaxis],
         ay[newaxis,:,newaxis,newaxis],
         bx[newaxis,newaxis,:,newaxis],
         by[newaxis,newaxis,newaxis,:])

or

p1 = np.array([[0,1], [1,3], [5,6], [7,9]])
p2 = np.array([[3,6], [2,7], [9,2], [1,3]])

def fun(a, b):
     return np.hypot(a[...,0] - b[...,0], a[...,1] - b[...,1])

d = fun(p1[:,newaxis,:], p2[newaxis,:,:])

depending on what it is that you actually want to evaluate.

--
Pauli Virtanen
_______________________________________________
SciPy-User mailing list
[hidden email]
https://mail.python.org/mailman/listinfo/scipy-user
Reply | Threaded
Open this post in threaded view
|

Re: Creating meshgrid from meshgrid

Florian Lindner
Hi,

thanks for your input!

Am 06.09.2017 um 16:05 schrieb Pauli Virtanen:

> Florian Lindner kirjoitti 06.09.2017 klo 06:31:
> [clip]
>>      def fun(x):
>>          return x
>>
>>      mgrid = np.meshgrid( a,b )
>>      A = fun( np.abs(mgrid[0] - mgrid[1] ) )
>>
>> Result A is of size N x N:
>>
>>      array([[ 9,  8,  7],
>>             [19, 18, 17],
>>             [29, 28, 27]])
>
> Don't use meshgrid, use broadcasting:
>
> from numpy import newaxis
>
> ax = np.array([0, 1, 2])
> ay = np.array([4, 5])
> bx = np.array([7, 8, 9, 10])
> by = np.array([11, 12, 13, 14])
>
> def fun(x1, y1, x2, y2):
>     return abs(x1 - x2) + abs(y1 - y2)
>
> d = fun(ax[:,newaxis,newaxis,newaxis],
>         ay[newaxis,:,newaxis,newaxis],
>         bx[newaxis,newaxis,:,newaxis],
>         by[newaxis,newaxis,newaxis,:])

Ok, this gives me what I want after a d.reshape(6, 16).
However, it is very tied to two dimensions and seem to not easily generalizable to more (or less) dimensions.

> or
>
> p1 = np.array([[0,1], [1,3], [5,6], [7,9]])
> p2 = np.array([[3,6], [2,7], [9,2], [1,3]])
>
> def fun(a, b):
>     return np.hypot(a[...,0] - b[...,0], a[...,1] - b[...,1])
>
> d = fun(p1[:,newaxis,:], p2[newaxis,:,:])
>
> depending on what it is that you actually want to evaluate.

That returns d of shape (4, 4), where I expect (16, 16).

I try to reformulate my problem:

From:

ax = np.array([1,2,3])
ay = np.array([3,4,5])
a = np.meshgrid(ax, ay)

bx = np.array([10,20,30])
by = np.array([30,40,50])
b = np.meshgrid(bx, by)

I can build the coordinates of all points in the meshes:

ca = np.array(list(zip(a[0].flatten(), a[1].flatten())))
cb = np.array(list(zip(b[0].flatten(), b[1].flatten())))

(in 1-d it will be just like ca = np.linspace(...))

just like your p1 and p2.

Now:

A = np.zeros([len(ca), len(cb)])
for i, x in enumerate(ca):
    for j, y in enumerate(cb):
        A[i, j] = func(np.linalg.norm(x-y))

gives me exactly what I want and should work for any number of dimensions in ca/cb. , albeit in a very un-numpythonic
and probably slow way.

You have any hints on how to improve that?

Best Thanks,
Florian

_______________________________________________
SciPy-User mailing list
[hidden email]
https://mail.python.org/mailman/listinfo/scipy-user