An example from Plenoptic Modeling:An Image-Based Rendering System Eq. 22
The original equation:
I❤️LA implementation:
r̄ = v̄×ō
s̄ = ō×ū
n̄ = ū×v̄
`kᵣ` = r̄⋅(`C̄ₐ`-V̄)
`kₛ` = s̄⋅(`C̄ₐ`-V̄)
`kₙ` = n̄⋅(`C̄ₐ`-V̄)
`x(θ,v)` = (r̄⋅`D_A`(θ, v)+`kᵣ`δ(θ, v))/(n̄⋅`D_A`(θ, v)+`kₙ`δ(θ, v))
`y(θ,v)` = (s̄⋅`D_A`(θ, v)+`kₛ`δ(θ, v))/(n̄⋅`D_A`(θ, v)+`kₙ`δ(θ, v))
where
v̄ ∈ ℝ^3
ō ∈ ℝ^3
ū ∈ ℝ^3
V̄ ∈ ℝ^3
`C̄ₐ` ∈ ℝ^3
θ ∈ ℝ
v ∈ ℝ
`D_A`: ℝ,ℝ → ℝ^3
δ: ℝ,ℝ → ℝ
I❤️LA compiled to C++/Eigen:
/*
r̄ = v̄×ō
s̄ = ō×ū
n̄ = ū×v̄
`kᵣ` = r̄⋅(`C̄ₐ`-V̄)
`kₛ` = s̄⋅(`C̄ₐ`-V̄)
`kₙ` = n̄⋅(`C̄ₐ`-V̄)
`x(θ,v)` = (r̄⋅`D_A`(θ, v)+`kᵣ`δ(θ, v))/(n̄⋅`D_A`(θ, v)+`kₙ`δ(θ, v))
`y(θ,v)` = (s̄⋅`D_A`(θ, v)+`kₛ`δ(θ, v))/(n̄⋅`D_A`(θ, v)+`kₙ`δ(θ, v))
where
v̄ ∈ ℝ^3
ō ∈ ℝ^3
ū ∈ ℝ^3
V̄ ∈ ℝ^3
`C̄ₐ` ∈ ℝ^3
θ ∈ ℝ
v ∈ ℝ
`D_A`: ℝ,ℝ → ℝ^3
δ: ℝ,ℝ → ℝ
*/
#include <Eigen/Core>
#include <Eigen/Dense>
#include <Eigen/Sparse>
#include <iostream>
#include <set>
struct plenoptic_modeling_22ResultType {
Eigen::Matrix<double, 3, 1> r̄;
Eigen::Matrix<double, 3, 1> s̄;
Eigen::Matrix<double, 3, 1> n̄;
double kᵣ;
double kₛ;
double kₙ;
double x_left_parenthesis_θ_comma_v_right_parenthesis;
double y_left_parenthesis_θ_comma_v_right_parenthesis;
plenoptic_modeling_22ResultType(const Eigen::Matrix<double, 3, 1> & r̄,
const Eigen::Matrix<double, 3, 1> & s̄,
const Eigen::Matrix<double, 3, 1> & n̄,
const double & kᵣ,
const double & kₛ,
const double & kₙ,
const double & x_left_parenthesis_θ_comma_v_right_parenthesis,
const double & y_left_parenthesis_θ_comma_v_right_parenthesis)
: r̄(r̄),
s̄(s̄),
n̄(n̄),
kᵣ(kᵣ),
kₛ(kₛ),
kₙ(kₙ),
x_left_parenthesis_θ_comma_v_right_parenthesis(x_left_parenthesis_θ_comma_v_right_parenthesis),
y_left_parenthesis_θ_comma_v_right_parenthesis(y_left_parenthesis_θ_comma_v_right_parenthesis)
{}
};
/**
* plenoptic_modeling_22
*
* @param D_A ℝ,ℝ → ℝ^3
* @param δ ℝ,ℝ → ℝ
* @return y_left_parenthesis_θ_comma_v_right_parenthesis
*/
plenoptic_modeling_22ResultType plenoptic_modeling_22(
const Eigen::Matrix<double, 3, 1> & v̄,
const Eigen::Matrix<double, 3, 1> & ō,
const Eigen::Matrix<double, 3, 1> & ū,
const Eigen::Matrix<double, 3, 1> & V̄,
const Eigen::Matrix<double, 3, 1> & C_combining_macron_ₐ,
const double & θ,
const double & v,
const std::function<Eigen::Matrix<double, 3, 1>(double, double)> & D_A,
const std::function<double(double, double)> & δ)
{
Eigen::Matrix<double, 3, 1> r̄ = (v̄).cross(ō);
Eigen::Matrix<double, 3, 1> s̄ = (ō).cross(ū);
Eigen::Matrix<double, 3, 1> n̄ = (ū).cross(v̄);
double kᵣ = (r̄).dot((C_combining_macron_ₐ - V̄));
double kₛ = (s̄).dot((C_combining_macron_ₐ - V̄));
double kₙ = (n̄).dot((C_combining_macron_ₐ - V̄));
double x_left_parenthesis_θ_comma_v_right_parenthesis = ((r̄).dot(D_A(θ, v)) + kᵣ * δ(θ, v)) / double(((n̄).dot(D_A(θ, v)) + kₙ * δ(θ, v)));
double y_left_parenthesis_θ_comma_v_right_parenthesis = ((s̄).dot(D_A(θ, v)) + kₛ * δ(θ, v)) / double(((n̄).dot(D_A(θ, v)) + kₙ * δ(θ, v)));
return plenoptic_modeling_22ResultType(r̄, s̄, n̄, kᵣ, kₛ, kₙ, x_left_parenthesis_θ_comma_v_right_parenthesis, y_left_parenthesis_θ_comma_v_right_parenthesis);
}
void generateRandomData(Eigen::Matrix<double, 3, 1> & v̄,
Eigen::Matrix<double, 3, 1> & ō,
Eigen::Matrix<double, 3, 1> & ū,
Eigen::Matrix<double, 3, 1> & V̄,
Eigen::Matrix<double, 3, 1> & C_combining_macron_ₐ,
double & θ,
double & v,
std::function<Eigen::Matrix<double, 3, 1>(double, double)> & D_A,
std::function<double(double, double)> & δ)
{
θ = rand() % 10;
v = rand() % 10;
v̄ = Eigen::VectorXd::Random(3);
ō = Eigen::VectorXd::Random(3);
ū = Eigen::VectorXd::Random(3);
V̄ = Eigen::VectorXd::Random(3);
C_combining_macron_ₐ = Eigen::VectorXd::Random(3);
D_A = [](double, double)->Eigen::Matrix<double, 3, 1>{
return Eigen::VectorXd::Random(3);
};
δ = [](double, double)->double{
return rand() % 10;
};
}
int main(int argc, char *argv[])
{
srand((int)time(NULL));
Eigen::Matrix<double, 3, 1> v̄;
Eigen::Matrix<double, 3, 1> ō;
Eigen::Matrix<double, 3, 1> ū;
Eigen::Matrix<double, 3, 1> V̄;
Eigen::Matrix<double, 3, 1> C_combining_macron_ₐ;
double θ;
double v;
std::function<Eigen::Matrix<double, 3, 1>(double, double)> D_A;
std::function<double(double, double)> δ;
generateRandomData(v̄, ō, ū, V̄, C_combining_macron_ₐ, θ, v, D_A, δ);
plenoptic_modeling_22ResultType func_value = plenoptic_modeling_22(v̄, ō, ū, V̄, C_combining_macron_ₐ, θ, v, D_A, δ);
std::cout<<"return value:\n"<<func_value.y_left_parenthesis_θ_comma_v_right_parenthesis<<std::endl;
return 0;
}
I❤️LA compiled to Python/NumPy/SciPy:
"""
r̄ = v̄×ō
s̄ = ō×ū
n̄ = ū×v̄
`kᵣ` = r̄⋅(`C̄ₐ`-V̄)
`kₛ` = s̄⋅(`C̄ₐ`-V̄)
`kₙ` = n̄⋅(`C̄ₐ`-V̄)
`x(θ,v)` = (r̄⋅`D_A`(θ, v)+`kᵣ`δ(θ, v))/(n̄⋅`D_A`(θ, v)+`kₙ`δ(θ, v))
`y(θ,v)` = (s̄⋅`D_A`(θ, v)+`kₛ`δ(θ, v))/(n̄⋅`D_A`(θ, v)+`kₙ`δ(θ, v))
where
v̄ ∈ ℝ^3
ō ∈ ℝ^3
ū ∈ ℝ^3
V̄ ∈ ℝ^3
`C̄ₐ` ∈ ℝ^3
θ ∈ ℝ
v ∈ ℝ
`D_A`: ℝ,ℝ → ℝ^3
δ: ℝ,ℝ → ℝ
"""
import numpy as np
import scipy
import scipy.linalg
from scipy import sparse
from scipy.integrate import quad
from scipy.optimize import minimize
class plenoptic_modeling_22ResultType:
def __init__( self, r̄, s̄, n̄, kᵣ, kₛ, kₙ, x_left_parenthesis_θ_comma_v_right_parenthesis, y_left_parenthesis_θ_comma_v_right_parenthesis):
self.r̄ = r̄
self.s̄ = s̄
self.n̄ = n̄
self.kᵣ = kᵣ
self.kₛ = kₛ
self.kₙ = kₙ
self.x_left_parenthesis_θ_comma_v_right_parenthesis = x_left_parenthesis_θ_comma_v_right_parenthesis
self.y_left_parenthesis_θ_comma_v_right_parenthesis = y_left_parenthesis_θ_comma_v_right_parenthesis
def plenoptic_modeling_22(v̄, ō, ū, V̄, C_combining_macron_ₐ, θ, v, D_A, δ):
"""
:param :D_A : ℝ,ℝ → ℝ^3
:param :δ : ℝ,ℝ → ℝ
"""
v̄ = np.asarray(v̄, dtype=np.float64)
ō = np.asarray(ō, dtype=np.float64)
ū = np.asarray(ū, dtype=np.float64)
V̄ = np.asarray(V̄, dtype=np.float64)
C_combining_macron_ₐ = np.asarray(C_combining_macron_ₐ, dtype=np.float64)
assert v̄.shape == (3,)
assert ō.shape == (3,)
assert ū.shape == (3,)
assert V̄.shape == (3,)
assert C_combining_macron_ₐ.shape == (3,)
assert np.ndim(θ) == 0
assert np.ndim(v) == 0
r̄ = np.cross(v̄, ō)
s̄ = np.cross(ō, ū)
n̄ = np.cross(ū, v̄)
kᵣ = np.dot((r̄).ravel(), ((C_combining_macron_ₐ - V̄)).ravel())
kₛ = np.dot((s̄).ravel(), ((C_combining_macron_ₐ - V̄)).ravel())
kₙ = np.dot((n̄).ravel(), ((C_combining_macron_ₐ - V̄)).ravel())
x_left_parenthesis_θ_comma_v_right_parenthesis = (np.dot((r̄).ravel(), (D_A(θ, v)).ravel()) + kᵣ * δ(θ, v)) / (np.dot((n̄).ravel(), (D_A(θ, v)).ravel()) + kₙ * δ(θ, v))
y_left_parenthesis_θ_comma_v_right_parenthesis = (np.dot((s̄).ravel(), (D_A(θ, v)).ravel()) + kₛ * δ(θ, v)) / (np.dot((n̄).ravel(), (D_A(θ, v)).ravel()) + kₙ * δ(θ, v))
return plenoptic_modeling_22ResultType(r̄, s̄, n̄, kᵣ, kₛ, kₙ, x_left_parenthesis_θ_comma_v_right_parenthesis, y_left_parenthesis_θ_comma_v_right_parenthesis)
def generateRandomData():
θ = np.random.randn()
v = np.random.randn()
v̄ = np.random.randn(3)
ō = np.random.randn(3)
ū = np.random.randn(3)
V̄ = np.random.randn(3)
C_combining_macron_ₐ = np.random.randn(3)
def D_A(p0, p1):
return np.random.randn(3)
def δ(p0, p1):
return np.random.randn()
return v̄, ō, ū, V̄, C_combining_macron_ₐ, θ, v, D_A, δ
if __name__ == '__main__':
v̄, ō, ū, V̄, C_combining_macron_ₐ, θ, v, D_A, δ = generateRandomData()
print("v̄:", v̄)
print("ō:", ō)
print("ū:", ū)
print("V̄:", V̄)
print("C_combining_macron_ₐ:", C_combining_macron_ₐ)
print("θ:", θ)
print("v:", v)
print("D_A:", D_A)
print("δ:", δ)
func_value = plenoptic_modeling_22(v̄, ō, ū, V̄, C_combining_macron_ₐ, θ, v, D_A, δ)
print("return value: ", func_value.y_left_parenthesis_θ_comma_v_right_parenthesis)
I❤️LA compiled to MATLAB:
function output = plenoptic_modeling_22(v_bar, o_bar, u_bar, V_bar, C_bar_a, theta, v, D_A, delta)
% output = plenoptic_modeling_22(v̄, ō, ū, V̄, `C̄ₐ`, θ, v, `D_A`, δ)
%
% r̄ = v̄×ō
% s̄ = ō×ū
% n̄ = ū×v̄
%
% `kᵣ` = r̄⋅(`C̄ₐ`-V̄)
% `kₛ` = s̄⋅(`C̄ₐ`-V̄)
% `kₙ` = n̄⋅(`C̄ₐ`-V̄)
%
% `x(θ,v)` = (r̄⋅`D_A`(θ, v)+`kᵣ`δ(θ, v))/(n̄⋅`D_A`(θ, v)+`kₙ`δ(θ, v))
% `y(θ,v)` = (s̄⋅`D_A`(θ, v)+`kₛ`δ(θ, v))/(n̄⋅`D_A`(θ, v)+`kₙ`δ(θ, v))
%
% where
%
% v̄ ∈ ℝ^3
% ō ∈ ℝ^3
% ū ∈ ℝ^3
% V̄ ∈ ℝ^3
% `C̄ₐ` ∈ ℝ^3
% θ ∈ ℝ
% v ∈ ℝ
% `D_A`: ℝ,ℝ → ℝ^3
% δ: ℝ,ℝ → ℝ
if nargin==0
warning('generating random input data');
[v_bar, o_bar, u_bar, V_bar, C_bar_a, theta, v, D_A, delta] = generateRandomData();
end
function [v_bar, o_bar, u_bar, V_bar, C_bar_a, theta, v, D_A, delta] = generateRandomData()
theta = randn();
v = randn();
v_bar = randn(3,1);
o_bar = randn(3,1);
u_bar = randn(3,1);
V_bar = randn(3,1);
C_bar_a = randn(3,1);
D_A = @D_AFunc;
rseed = randi(2^32);
function tmp = D_AFunc(p0, p1)
rng(rseed);
tmp = randn(3,1);
end
delta = @deltaFunc;
rseed = randi(2^32);
function tmp = deltaFunc(p0, p1)
rng(rseed);
tmp = randn();
end
end
v_bar = reshape(v_bar,[],1);
o_bar = reshape(o_bar,[],1);
u_bar = reshape(u_bar,[],1);
V_bar = reshape(V_bar,[],1);
C_bar_a = reshape(C_bar_a,[],1);
assert( numel(v_bar) == 3 );
assert( numel(o_bar) == 3 );
assert( numel(u_bar) == 3 );
assert( numel(V_bar) == 3 );
assert( numel(C_bar_a) == 3 );
assert(numel(theta) == 1);
assert(numel(v) == 1);
r_bar = cross(v_bar, o_bar);
s_bar = cross(o_bar, u_bar);
n_bar = cross(u_bar, v_bar);
k_r = dot(r_bar,(C_bar_a - V_bar));
k_s = dot(s_bar,(C_bar_a - V_bar));
k_n = dot(n_bar,(C_bar_a - V_bar));
x_theta_v = (dot(r_bar,D_A(theta, v)) + k_r * delta(theta, v)) / (dot(n_bar,D_A(theta, v)) + k_n * delta(theta, v));
y_theta_v = (dot(s_bar,D_A(theta, v)) + k_s * delta(theta, v)) / (dot(n_bar,D_A(theta, v)) + k_n * delta(theta, v));
output.r_bar = r_bar;
output.s_bar = s_bar;
output.n_bar = n_bar;
output.k_r = k_r;
output.k_s = k_s;
output.k_n = k_n;
output.x_theta_v = x_theta_v;
output.y_theta_v = y_theta_v;
end
I❤️LA compiled to LaTeX:
\documentclass[12pt]{article}
\usepackage{mathdots}
\usepackage[bb=boondox]{mathalfa}
\usepackage{mathtools}
\usepackage{amssymb}
\usepackage{libertine}
\DeclareMathOperator*{\argmax}{arg\,max}
\DeclareMathOperator*{\argmin}{arg\,min}
\usepackage[paperheight=8in,paperwidth=4in,margin=.3in,heightrounded]{geometry}
\let\originalleft\left
\let\originalright\right
\renewcommand{\left}{\mathopen{}\mathclose\bgroup\originalleft}
\renewcommand{\right}{\aftergroup\egroup\originalright}
\begin{document}
\begin{center}
\resizebox{\textwidth}{!}
{
\begin{minipage}[c]{\textwidth}
\begin{align*}
\textit{r̄} & = \textit{v̄} × \textit{ō} \\
\textit{s̄} & = \textit{ō} × \textit{ū} \\
\textit{n̄} & = \textit{ū} × \textit{v̄} \\
\textit{k\textsubscript{r}} & = \textit{r̄} \cdot \left( \textit{C̄ₐ} - \textit{V̄} \right) \\
\textit{kₛ} & = \textit{s̄} \cdot \left( \textit{C̄ₐ} - \textit{V̄} \right) \\
\textit{kₙ} & = \textit{n̄} \cdot \left( \textit{C̄ₐ} - \textit{V̄} \right) \\
\textit{x(θ,v)} & = \frac{\textit{r̄} \cdot \textit{D\_A}\left( \mathit{θ},\mathit{v} \right) + \textit{k\textsubscript{r}}\mathit{δ}\left( \mathit{θ},\mathit{v} \right)}{\textit{n̄} \cdot \textit{D\_A}\left( \mathit{θ},\mathit{v} \right) + \textit{kₙ}\mathit{δ}\left( \mathit{θ},\mathit{v} \right)} \\
\textit{y(θ,v)} & = \frac{\textit{s̄} \cdot \textit{D\_A}\left( \mathit{θ},\mathit{v} \right) + \textit{kₛ}\mathit{δ}\left( \mathit{θ},\mathit{v} \right)}{\textit{n̄} \cdot \textit{D\_A}\left( \mathit{θ},\mathit{v} \right) + \textit{kₙ}\mathit{δ}\left( \mathit{θ},\mathit{v} \right)} \\
\intertext{where}
\textit{v̄} & \in \mathbb{R}^{ 3} \\
\textit{ō} & \in \mathbb{R}^{ 3} \\
\textit{ū} & \in \mathbb{R}^{ 3} \\
\textit{V̄} & \in \mathbb{R}^{ 3} \\
\textit{C̄ₐ} & \in \mathbb{R}^{ 3} \\
\mathit{θ} & \in \mathbb{R} \\
\mathit{v} & \in \mathbb{R} \\
\textit{D\_A} & \in \mathbb{R},\mathbb{R}\rightarrow \mathbb{R}^{ 3} \\
\mathit{δ} & \in \mathbb{R},\mathbb{R}\rightarrow \mathbb{R} \\
\\
\end{align*}
\end{minipage}
}
\end{center}
\end{document}
I❤️LA LaTeX output: