How do I numerically calculate a function from its gradient? I know the gradient of a function t on a cartesian grid:
$\vec{g}(xi,yj,zk)=\nabla t(xi,yj,zk)$.
I know t for the center pillar: $\ t(xc,yc,zk)$.
For each node in the cartesian grid I want to calculate the timeshift:
$\ ts(xi,yj,zk) = t(xi,yj,zk) - t(xc,yc,zk)$  
According to the gradient theorem: 
$\ ts(xi,yj,zk) = \int_{(xc,yc,zk)}^{(xi,yj,zk)} \vec{g}(\vec{r})d\vec{r}$
But how do I quck and runtime effective code this (in Matlab)?
I want to avoid summing over n contributions for every node in the grid; O(n^4) :-(. 
Thanks in advance for any answers! 
 A: You can do it in linear time by marching outwards along the grid starting from $(x_c,y_c,z_k)$ and filling in the unknown values one at a time. For example, at the node immediately to the right of the pillar, you have
$$\begin{align}
t(x_c+1,y_c,z_k) - t(x_c,y_c,z_k) &= \int_{(x_c,y_c,z_k)}^{(x_c+1,y_c,z_k)}\nabla t(\vec r)\cdot\mathrm d\vec r \\
&\approx \frac{g_x(x_c,y_c,z_k) + g_x(x_c+1,y_c,z_k)}2
\end{align}$$
and since you know all the other terms in the equation, you can fill in the value of $t(x_c+1,y_c,z_k)$. Now knowing $t(x_c+1,y_c,z_k)$, you can fill in the ones next to it, i.e. $t(x_c+2,y_c,z_k)$, $t(x_c+1,y_c+1,z_k)$, and so on, and continue similarly from there.
A: How about using the Fourier transform ?:
$\ [G_x, G_y, G_z] = i*[k_x*T_x, k_y*T_y, k_z*T_z] $
where $\ G_x $ is the x component of the Fourier transform of g and $\ T_x $ is the x component of the Fourier transform of t.
Solving for T:
$\ [T_x, T_y, T_z] = 1/i*[G_x/k_x, G_y/k_y, G_z/k_z] $
Find the kx, ky, kz from Nyquist: $\ \Delta k_x = \frac{1}{n_x*\Delta x} $, where $\ n_x $: the number of samples in x direction and $\ \Delta x $: the sampling in x direction.
To get the correct DC offset I am guessing that kx = ky = 0 should be placed in the position of xc, yc? 
Finally:
$\ t(x,y,z) = IFFT(T) $
Will this work?
