MATLAB собирается из корзины для Pol обратно в корзину для цилиндрических координат участка

Question

1. То, что я хочу сделать:
(я) Используйте вход п для создания п * п декартовой сетки

[x y] = meshgrid(linspace(-1,1,n)); 

(II) Генерировать полярные координаты

[theta r] = cart2pol(x,y); 

(III) Оценка фу nction в цилиндрических координатах

z = f(theta,r); 

(IV) Plot результат, используя (скажем) pcolor (или серфинг, или что-нибудь)

pcolor(x,y,abs(z).^2) %The function is complex, a Laguerre-Gauss to be exact. 

2. Что я могу сделать ... Только как я могу получить участки для работы, начиная с моими полярными параметрами и работать обратно в декартовую оттуда:
(я) Определения параметров

r=linspace(0,1,n); theta=linspace(0,2*pi,n); 

(б) Создать обе сетки и оценить е

[theta r]=meshgrid(theta,r); 
[x y]=pol2cart(theta,r); 
z=f(theta,r); 

(III) Участок

pcolor(x,y,abs(z).^2) 

Проблема заключается в том, что теперь моя решетка круговая, и я хотел бы, чтобы оценить функция всюду ВКЛЮЧЕНА ПРЯМОУГОЛЬНОЙ сеткой (потому что мой анализ зависит от наличия квадратных пиксельных массивов). Повторяя, используя метод 2 выше, я получаю круговой график, описанный в квадрате; представьте себе черный круг с белым по краям ... но я ХОЧУ оценить функцию в этой «белой» области. ОДНАКО, использование метода 1 НЕ работает - функция все испортилась, когда я завожу (просто режимы Google Laguerre-Gauss, чтобы посмотреть, как будут выглядеть сюжеты).

Я хочу, чтобы начать с прямоугольной сетки и назначить каждую точку полярной координаты, вместо того, чтобы начинать с полярных координат и назначать им все декартовы точки.

Я уже давно занимаюсь этим, и я не могу понять, как обойти эту, казалось бы, простую проблему.

Edit 1
кажется, что проблема заключается в том, как создаются координатные матрицы. Ниже я опубликовал скриншоты простого примера 3by3, иллюстрирующего, как подход 1 и подход 2 генерируют разные числа.

Как сделать эти номера совместимыми?

У меня нет точек репутации, поэтому я не могу загрузить изображения напрямую ... ссылки ниже показывают пример 3by3 ... см. Комментарии для ссылок на фактические изображения участков Laguerre-Gauss, которые я пытаюсь сделать ...

apply cart2pol
apply pol2cart

Edit 2

В настоящее время, в результате подхода (1.) дает неверные результаты, как показано здесь:

desired approach, wrong result

второй подход дает r но, к сожалению, это всего лишь круг, а не весь квадрат. Показано здесь:

implemented approach, limited result

3D графики обоих подходов приведены here - только красочная часть верхней фигуры правильно.

Редактировать 3

Вот screenshot of the function f, который используется выше. Обратите внимание, что он запрашивает больше входных параметров, чем только r,theta. Типичные значения:

w0 = 0.5; 
p = 0; 
l = 5; 

Функция C дает нормализацию и L являются полиномы Лагерра. Обе эти функции были тщательно протестированы и дали ожидаемые результаты.

Edit 4
Здесь достаточно кода, чтобы запустить мой пример z=U(0,5,r,phi,w0)+U(0,-5,r,phi,w0); явно. Сам сюжет дается pcolor(x,y,abs(z).^2).

Обратите внимание, что в качестве комментария добавлена ​​функция Lpl(). Это нужно будет сохранить как собственный m-файл для правильной работы U-функции.

%% Laguerre-Gauss Modes U = U(p,l,r,phi,w0) 
% Source: OAM theory paper section 2.A eqn 1. 
% Assuming POLAR coordinates and evaluating AT beam waist. 
% -- That is, z=0 for w(z)=w0(sqrt(1+z/zR)) 
% ---- ie, w(0) = w0 
% Assuming z=0 also renders the Gouy phase arctan(z/zR) irrelevant. 
% Note: Rayleigh Range zR is not explicitly defined because z=0 --> it is irrelevant too. 
% Since zR is the only wavelength dependent term, wavelength also doesn't 
% matter. 

function out = U(p,l,r,phi,w0) 
%Function handles for clarity 
e = @(x) exp(x); 
C = @(p,l) sqrt((2*factorial(p))/(pi*factorial(p+abs(l)))); 
L = @(p,l,z) Lpl(p,l,z); 

%% Lpl() FUNCTION 
% function out = Lpl(p,l,z) 
% 
% l=abs(l); 
% LL=0; 
% for mm=1:p+1 
%  m=mm-1; 
%  L=LL; 
%  LL= L+((-1)^m)*(factorial(p+l)/(factorial(p-m)*factorial(l+m)*factorial(m)))*(z.^m); 
% end 
% out = LL; 

%% 

out = (C(p,l)/w0)*... 
    (((sqrt(2).*r)/w0)^abs(l))*... 
    (e((-r.^2)/w0^2))*... 
    (L(p,l,((2.*r.^2)/w0^2)))*... 
    (e((-1)*1i*l.*phi)); `` 

Answer 1

Редактировать
Ответ был переписан на основе кода, представленной в Edit 4 вопроса.

Я думаю, что проблема связана с функцией U. Вы не применяете element wise operations ко всем частям уравнения. Если вы измените его на:

out = (C(p,l)./w0).* ...    % here it's a .* instead of * 
    (((sqrt(2).*r)./w0).^abs(l)).* ... % here it's a .* instead of * 
    (e((-r.^2)./w0.^2)).* ...   % here it's a .* instead of * 
    (L(p,l,((2.*r.^2)./w0.^2))).* ... % here it's a .* instead of * 
    (e((-1)*1i*l.*phi)); 

Вы получаете следующие два результата, показанные ниже.

Этот показатель используется ввод декартовых координат:

И эта цифра используется полярные координаты:

«грубее» Разрешение на втором рисунке обусловлено к менее подходящему разрешению сетки.Но по сути вы решаете те же функции.