二维legendre多项式

Legendre 多项式常用来表征方形波前的畸变。

一维legendre多项式

Legendre 多项式是一组在区间 $[-1,1]$ 上定义的多项式序列，它们具有正交性。具体来说，如果$P_n(x)$表示第 $n$阶 Legendre 多项式有：

正交性

对于任何两个不同的整数 $m和n$，有：

$${\int }_{-1}^1{P}_m(x)P_n(x)dx=0$$

这意味着不同阶次的 Legendre 多项式在这个区间上关于标准内积是相互正交的。

自正交性

对于同一个 $n$值，Legendre 多项式的自正交条件为：

$${\int }_{-1}^1[P_n(x){]}^{2}dx=\frac{2}{2n+1}$$

这些性质使得 Legendre 多项式在数学物理、工程学以及数值分析等领域中有着广泛的应用，特别是在解决涉及球对称性的物理问题时，如电磁场理论、量子力学中的角动量理论等。

二维Legendre多项式

二维Legendre多项式通常不是以单一形式定义的，而是通过组合一维Legendre多项式来构造。在二维空间中，我们可以通过将一维Legendre多项式 $P_n(x)$和$P_m(y)$组合起来，形成一个二维多项式 $P_n(x)P_m(y)$，这样可以保持正交性。

考虑在单位方形 [-1, 1] ×[-1, 1]上定义的二维Legendre多项式$P_{nm}(x,y)=P_n(x)P_m(y)$，其中 $P_n(x)$ 和 $P_m(y)$ 分别是 $x$ 和 $y$方向上的一维Legendre多项式。

正交性证明

为了证明$P_{nm}(x,y)$的正交性，我们需要计算两个不同阶次的二维Legendre多项式之间的内积：

$${\iint }_{[-1,1{]}^2}{P}_{nm}(x,y)P_{kl}(x,y)dxdy$$

这里 $P_{nm}(x,y)=P_n(x)P_m(y)$和$P_{kl}(x,y)=P_k(x)P_l(y)$将积分展开：

$${\iint }_{[-1,1{]}^2}{P}_n(x)P_m(y)P_k(x)P_l(y)dxdy=\left({\int }_{-1}^1{P}_n(x)P_k(x)dx\right)\left({\int }_{-1}^1{P}_m(y)P_l(y)dy\right)$$

根据一维Legendre多项式的正交性：

• 如果$n\ne k$，则 ${\int }_{-1}^1{P}_n(x)P_k(x)dx=0$；
• 如果$m\ne l$，则 ${\int }_{-1}^1{P}_m(y)P_l(y)dy=0$。

因此，当 $n\ne k$ 或$m\ne l$ 时，上述二重积分的结果为零，这表明 $P_{nm}(x,y)$ 和$P_{kl}(x,y)$是正交的。

当 $n=k$且 $m=l$时，即两个多项式相同的情况下，积分结果为：

$$\left({\int }_{-1}^1[P_n(x){]}^{2}dx\right)\left({\int }_{-1}^1[P_m(y){]}^{2}dy\right)=\left(\frac{2}{2n+1}\right)\left(\frac{2}{2m+1}\right)=\frac{4}{(2n+1)(2m+1)}$$

这就证明了二维Legendre多项式$P_{nm}(x,y)$在单位方形$[-1,1]\times [-1,1]$上是正交的，并且当 $n=k$和$m=l$时，它们的内积为 $\frac{4}{(2n+1)(2m+1)}$。这种构造方法保证了二维Legendre多项式在多变量函数的正交分解中具有重要的应用价值。

可视化二维 Legendre 多项式

在 MATLAB 中可视化二维 Legendre 多项式可以通过以下步骤实现。我们将使用 meshgrid 函数生成网格点，然后计算每个网格点上的多项式值，最后使用 surf 或 contour 函数进行可视化。

下面是一个具体的例子，假设我们要可视化 $P_2(x)P_3(y)$这个二维 Legendre 多项式：

1. 定义一维 Legendre 多项式：
   我们可以使用递归公式或内置函数来定义一维 Legendre 多项式。

2. 生成网格点：
   使用 meshgrid 函数生成$x$ 和 $y$ 的网格点。

3. 计算多项式值：
   在每个网格点上计算$P_2(x)$和 $P_3(y)$的值，然后乘积得到 $P_2(x)P_3(y)$。

4. 可视化：
   使用 surf 或 contour 函数绘制三维表面图或等高线图。

以下是完整的 MATLAB 代码：

% 定义一维 Legendre 多项式函数
function P = legendre_poly(n, x)
    if n == 0
        P = ones(size(x));
    elseif n == 1
        P = x;
    else
        P0 = ones(size(x));
        P1 = x;
        for k = 2:n
            P = ((2*k-1)*x.*P1 - (k-1)*P0) / k;
            P0 = P1;
            P1 = P;
        end
    end
end

% 生成网格点
[x, y] = meshgrid(linspace(-1, 1, 100), linspace(-1, 1, 100));

% 计算 P2(x) 和 P3(y)
P2_x = legendre_poly(2, x);
P3_y = legendre_poly(3, y);

% 计算二维 Legendre 多项式 P2(x) * P3(y)
P23 = P2_x .* P3_y;

% 绘制三维表面图
figure;
surf(x, y, P23);
shading interp;
colorbar;
xlabel('x');
ylabel('y');
zlabel('P_2(x) * P_3(y)');
title('2D Legendre Polynomial P_2(x) * P_3(y)');
view(3); % 设置视角为三维

% 绘制等高线图
figure;
contourf(x, y, P23, 20); % 20 表示等高线的数量
colorbar;
xlabel('x');
ylabel('y');
title('Contour Plot of 2D Legendre Polynomial P_2(x) * P_3(y)');

解释

1. 定义一维 Legendre 多项式函数：legendre_poly 函数使用递归公式计算 (n) 阶 Legendre 多项式。
2. 生成网格点：meshgrid 函数生成 $x$ 和$y$的网格点。
3. 计算多项式值：分别计算 $P_2(x)$和$P_3(y)$ 的值，然后乘积得到$P_2(x)\cdot P_3(y)$。
4. 可视化：使用 surf 函数绘制三维表面图，使用 contourf 函数绘制等高线图。

图片出自：《基于单帧焦面图像的波前相位反演方法研究》

Legendre拟合

要计算波前$P(x,y)$的 Legendre 系数，可以利用二维 Legendre 多项式的正交性。假设波前$P(x,y)$在单位方形$[-1,1]\times [-1,1]$上定义，我们可以将其展开为二维 Legendre 多项式的级数形式：

$$P(x,y)=\sum _{n=0}^{\infty }\sum _{m=0}^{\infty }{c}_{nm}{P}_n(x)P_m(y)$$

其中$P_n(x)$和$P_m(y)$是一维 Legendre 多项式，$c_{nm}$是对应的 Legendre 系数。

方法1

利用正交性，可以通过以下步骤计算$c_{nm}$：

1. 定义一维 Legendre 多项式：
   使用递归公式或内置函数来定义一维 Legendre 多项式。

2. 计算内积：
   利用正交性，计算$c_{nm}$的公式为：

$$c_{nm}=\frac{(2n+1)(2m+1)}{4}{\int }_{-1}^1{\int }_{-1}^{1}P(x,y)P_n(x)P_m(y)dxdy$$

MATLAB 实现

以下是一个 MATLAB 代码示例，用于计算波前 ( P(x, y) ) 的 Legendre 系数：

% 定义一维 Legendre 多项式函数
function P = legendre_poly(n, x)
    if n == 0
        P = ones(size(x));
    elseif n == 1
        P = x;
    else
        P0 = ones(size(x));
        P1 = x;
        for k = 2:n
            P = ((2*k-1)*x.*P1 - (k-1)*P0) / k;
            P0 = P1;
            P1 = P;
        end
    end
end

% 定义波前 P(x, y)
P = @(x, y) sin(pi*x).*cos(pi*y); % 示例波前

% 计算 Legendre 系数
N = 5; % 最大阶数
c = zeros(N+1, N+1);

for n = 0:N
    for m = 0:N
        % 生成网格点
        [X, Y] = meshgrid(linspace(-1, 1, 100), linspace(-1, 1, 100));
        
        % 计算 P_n(x) 和 P_m(y)
        Pn_x = legendre_poly(n, X);
        Pm_y = legendre_poly(m, Y);
        
        % 计算内积
        integrand = P(X, Y) .* Pn_x .* Pm_y;
        c(n+1, m+1) = (2*n+1)*(2*m+1)/4 * integral2(@(x, y) P(x, y) .* legendre_poly(n, x) .* legendre_poly(m, y), -1, 1, -1, 1);
    end
end

% 显示系数
disp('Legendre 系数矩阵:');
disp(c);

解释

1. 定义一维 Legendre 多项式函数：legendre_poly 函数使用递归公式计算 $n$ 阶 Legendre 多项式。
2. 定义波前 $P(x,y)$：这里使用了一个示例波前 $P(x,y)=\sin (\pi x)\cos (\pi y)$。
3. 计算 Legendre 系数：
   • 生成网格点。
   • 计算每个网格点上的$P_n(x)$和$P_m(y)$。
   • 计算内积${\int }_{-1}^1{\int }_{-1}^{1}P(x,y)P_n(x)P_m(y)dxdy$。
   • 利用正交性公式计算$c_{nm}$。

方法2

也可以可以借鉴 Zernike 系数计算的方法来计算 Legendre 系数。 Zernike-Polynomials-MATLAB。
以下是一个 MATLAB 代码示例，用于计算波前$P(x,y)$的 Legendre 系数。

1. 定义一维 Legendre 多项式函数

首先，我们需要定义一个函数来计算一维 Legendre 多项式。

function P = legendre_poly(n, x)
    if n == 0
        P = ones(size(x));
    elseif n == 1
        P = x;
    else
        P0 = ones(size(x));
        P1 = x;
        for k = 2:n
            P = ((2*k-1)*x.*P1 - (k-1)*P0) / k;
            P0 = P1;
            P1 = P;
        end
    end
end

2. 生成二维 Legendre 多项式矩阵

接下来，我们需要一个函数来生成二维 Legendre 多项式矩阵。

function legMats = legendre_mats(im, indices)
    dim = size(im);
    x = linspace(-1, 1, dim(2));
    y = linspace(-1, 1, dim(1));
    [X, Y] = meshgrid(x, y);
    
    legMats = [];
    
    for i = 1:size(indices, 1)
        n = indices(i, 1);
        m = indices(i, 2);
        Pn_x = legendre_poly(n, X);
        Pm_y = legendre_poly(m, Y);
        legMats(:,:,i) = Pn_x .* Pm_y;
    end
end

3. 计算 Legendre 系数

最后，我们编写一个函数来计算波前$P(x,y)$的 Legendre 系数。

function [legMoments] = legendre_moments(P, indices)
    dim = size(P);
    NumCols = dim(2);
    NumRows = dim(1);
    
    legMats = legendre_mats(P, indices);
    
    leg_mats_reshaped = [];
    
    for i = 1:size(indices, 1)
        leg_mats_reshaped(:, i) = reshape(legMats(:, :, i), NumCols * NumRows, 1);
    end
    
    P(isnan(P)) = 0;
    leg_mats_reshaped(isnan(leg_mats_reshaped)) = 0;
    
    P_reshaped = reshape(P, NumCols * NumRows, 1);
    A = (leg_mats_reshaped.' * leg_mats_reshaped) \ (leg_mats_reshaped.' * P_reshaped);
    legMoments = A;
end

4. 示例用法

以下是一个示例，展示如何使用这些函数来计算波前$P(x,y)$的 Legendre 系数。

% 创建一些波前 P(x, y) 过一个 100x100 网格
x = linspace(-1, 1, 100);
y = linspace(-1, 1, 100);
[X, Y] = meshgrid(x, y);
P = sin(pi * X) .* cos(pi * Y);

% 指定希望使用的 Legendre 多项式的阶数
indices = [];
for n = 0:3
    for m = 0:3
        indices = [indices; n m];
    end
end

% 计算 Legendre 系数
legMoments = legendre_moments(P, indices);

% 显示系数及其对应的 (n, m) 索引
disp(['    Coeff     ', 'n        ', 'm']);
disp([legMoments indices]);

解释

1. 定义一维 Legendre 多项式函数：legendre_poly 函数使用递归公式计算 $n$阶 Legendre 多项式。
2. 生成二维 Legendre 多项式矩阵：legendre_mats 函数生成指定阶数的二维 Legendre 多项式矩阵。
3. 计算 Legendre 系数：legendre_moments 函数通过最小二乘法计算波前 $P(x,y)$的 Legendre 系数。
4. 示例用法：创建一个示例波前$P(x,y)$，指定希望使用的 Legendre 多项式的阶数，并计算相应的 Legendre 系数。