meshgrid a spirale in MATLAB

Question

Sto provando a produrre alcuni ologrammi generati al computer usando MATLAB. Ho usato maglia della griglia equidistanziati per inizializzare la griglia spaziale, e ho ottenuto immagine

Questo modello è sorta di quello che serve eccetto la regione centrale. La frangia dovrebbe essere nitida ma sfocata. Penso che potrebbe essere il problema della rete a maglie. Ho provato a generare una griglia in coordinate polari e la mappa in coordinate cartesiane utilizzando la funzione pol2cart di MATLAB. Sfortunatamente, non funziona altrettanto bene. Si potrebbe suggerire l'uso di griglie fini. Non funziona anche Penso che se riesco a generare una griglia a maglie a spirale, forse il problema è risolvibile. Inoltre, il numero delle braccia a spirale potrebbe, in generale, essere arbitrario, qualcuno potrebbe darmi un suggerimento su questo?

Ho allegato il codice (i miei progetti finali non sono esattamente uguali, ma ha un problema simile).

clc; clear all; close all; 
%% initialization 
tic 
lambda = 1.55e-6; 
k0 = 2*pi/lambda; 
c0 = 3e8; 
eta0 = 377; 
scale = 0.25e-6; 
NELEMENTS = 1600; 
GoldenRatio = (1+sqrt(5))/2; 
g = 2*pi*(1-1/GoldenRatio); 

pntsrc = zeros(NELEMENTS, 3); 
phisrc = zeros(NELEMENTS, 1); 
for idxe = 1:NELEMENTS 
    pntsrc(idxe, :) = scale*sqrt(idxe)*[cos(idxe*g), sin(idxe*g), 0]; 
    phisrc(idxe) = angle(-sin(idxe*g)+1i*cos(idxe*g)); 
end 
phisrc = 3*phisrc/2; % 3 arms (topological charge ell=3) 

%% post processing 
sigma = 1; 
polfilter = [0, 0, 1i*sigma; 0, 0, -1; -1i*sigma, 1, 0]; % cp filter 

xboundl = -100e-6; xboundu = 100e-6; 
yboundl = -100e-6; yboundu = 100e-6; 
xf = linspace(xboundl, xboundu, 100); 
yf = linspace(yboundl, yboundu, 100); 
zf = -400e-6; 
[pntobsx, pntobsy] = meshgrid(xf, yf); 
% how to generate a right mesh grid such that we can generate a decent result? 
pntobs = [pntobsx(:), pntobsy(:), zf*ones(size(pntobsx(:)))]; 
% arbitrary mesh may result in "wrong" results 

NPNTOBS = size(pntobs, 1); 
nxp = length(xf); 
nyp = length(yf); 

%% observation 
Eobs = zeros(NPNTOBS, 3); 

matlabpool open local 12 
parfor nobs = 1:NPNTOBS 
    rp = pntobs(nobs, :); 
    Erad = [0; 0; 0]; 
    for idx = 1:NELEMENTS 
    rs = pntsrc(idx, :); 
    p = exp(sigma*1i*2*phisrc(idx))*[1 -sigma*1i 0]/2; % simplified here 
    u = rp - rs; 
    r = sqrt(u(1)^2+u(2)^2+u(3)^2); %norm(u); 
    u = u/r; % unit vector 
    ut = [u(2)*p(3)-u(3)*p(2),... 
     u(3)*p(1)-u(1)*p(3), ... 
     u(1)*p(2)-u(2)*p(1)]; % cross product: u cross p 
    Erad = Erad + ... % u cross p cross u, do not use the built-in func 
     c0*k0^2/4/pi*exp(1i*k0*r)/r*eta0*... 
     [ut(2)*u(3)-ut(3)*u(2);... 
     ut(3)*u(1)-ut(1)*u(3); ... 
     ut(1)*u(2)-ut(2)*u(1)]; 
    end 
    Eobs(nobs, :) = Erad; % filter neglected here 
end 
matlabpool close 
Eobs = Eobs/max(max(sum(abs(Eobs), 2))); % normailized 

%% source, gaussian beam 
E0 = 1; 
w0 = 80e-6; 
theta = 0; % may be titled 
RotateX = [1, 0, 0; ... 
    0, cosd(theta), -sind(theta); ... 
    0, sind(theta), cosd(theta)]; 

Esrc = zeros(NPNTOBS, 3); 
for nobs = 1:NPNTOBS 
    rp = RotateX*[pntobs(nobs, 1:2).'; 0]; 
    z = rp(3); 
    r = sqrt(sum(abs(rp(1:2)).^2)); 
    zR = pi*w0^2/lambda; 
    wz = w0*sqrt(1+z^2/zR^2); 
    Rz = z^2+zR^2; 
    zetaz = atan(z/zR); 
    gaussian = E0*w0/wz*exp(-r^2/wz^2-1i*k0*z-1i*k0*0*r^2/Rz/2+1i*zetaz);% ... 
    Esrc(nobs, :) = (polfilter*gaussian*[1; -1i; 0]).'/sqrt(2)/2; 
end 
Esrc = [Esrc(:, 2), Esrc(:, 3), Esrc(:, 1)]; 
Esrc = Esrc/max(max(sum(abs(Esrc), 2))); % normailized 
toc 

%% visualization 
fringe = Eobs + Esrc; % I'll have a different formula in my code 
normEsrc = reshape(sum(abs(Esrc).^2, 2), [nyp nxp]); 
normEobs = reshape(sum(abs(Eobs).^2, 2), [nyp nxp]); 
normFringe = reshape(sum(abs(fringe).^2, 2), [nyp nxp]); 

close all; 
xf0 = linspace(xboundl, xboundu, 500); 
yf0 = linspace(yboundl, yboundu, 500); 
[xfi, yfi] = meshgrid(xf0, yf0); 
data = interp2(xf, yf, normFringe, xfi, yfi); 
figure; surf(xfi, yfi, data,'edgecolor','none'); 
% tri = delaunay(xfi, yfi); trisurf(tri, xfi, yfi, data, 'edgecolor','none'); 
xlim([xboundl, xboundu]) 
ylim([yboundl, yboundu]) 
% colorbar 
view(0,90) 
colormap(hot) 
axis equal 
axis off 
title('fringe thereo. ', ... 
    'fontsize', 18) 

Answer 1

Non ho letto il codice perché è troppo lungo per fare una cosa così semplice. Ho scritto la mia e qui è il risultato:

il codice è

%spiral.m 
function val = spiral(x,y) 

    r = sqrt(x*x + y*y); 
    a = atan2(y,x)*2+r;     

    x = r*cos(a); 
    y = r*sin(a); 

    val = exp(-x*x*y*y); 

    val = 1/(1+exp(-1000*(val))); 

endfunction 


%show.m 
n=300; 
l = 7; 
A = zeros(n); 

for i=1:n 
for j=1:n 
    A(i,j) = spiral(2*(i/n-0.5)*l,2*(j/n-0.5)*l); 
end 
end 


imshow(A) %don't know if imshow is in matlab. I used octave. 

la chiave per i sharpnes è la linea

val = 1/(1+exp(-1000*(val))); 

È logistic function. Il numero 1000 definisce quanto sarà nitida la tua immagine. Quindi abbassalo per un'immagine più sfocata o più alto per renderlo più nitido.

Spero che questo risponde alla tua domanda;)

Edit: E 'vero divertimento con cui giocare. Ecco un altro spirale:

function val = spiral(x,y) 

    s= 0.5; 

    r = sqrt(x*x + y*y); 
    a = atan2(y,x)*2+r*r*r;     

    x = r*cos(a); 
    y = r*sin(a); 

    val = 0; 
    if (abs(x)<s) 
    val = s-abs(x); 
    endif 
if(abs(y)<s) 
    val =max(s-abs(y),val); 
    endif 

    %val = 1/(1+exp(-1*(val))); 

endfunction 

Edit2: Fun, fun, fun! Qui le braccia non si assottigliano.

function val = spiral(x,y) 

    s= 0.1; 

    r = sqrt(x*x + y*y); 
    a = atan2(y,x)*2+r*r;     % h 

    x = r*cos(a); 
    y = r*sin(a); 

    val = 0; 
    s = s*exp(r); 
    if (abs(x)<s) 
    val = s-abs(x); 
    endif 
if(abs(y)<s) 
    val =max(s-abs(y),val); 
    endif 
val = val/s; 
val = 1/(1+exp(-10*(val))); 

endfunction 

Accidenti tua domanda ho davvero bisogno di studiare per il mio esame, arghhh!

Edit3:

ho Vectorised il codice e corre molto più veloce.

%spiral.m 
    function val = spiral(x,y) 
    s= 2; 

    r = sqrt(x.*x + y.*y); 
    a = atan2(y,x)*8+exp(r);     

    x = r.*cos(a); 
    y = r.*sin(a); 

    val = 0; 
    s = s.*exp(-0.1*r); 
    val = r; 
    val = (abs(x)<s).*(s-abs(x)); 
    val = val./s; 

% val = 1./(1.+exp(-1*(val))); 
endfunction 


%show.m 

n=1000; 
l = 3; 
A = zeros(n); 
[X,Y] = meshgrid(-l:2*l/n:l); 
A = spiral(X,Y); 
imshow(A) 

Answer 2

Siamo spiacenti, non è possibile inserire cifre. Ma questo potrebbe aiutare. L'ho scritto per esperimenti con modulatori spaziali di ampiezza ...

R=70;   % radius of curvature of fresnel lens (in pixel units) 
A=0;    % oblique incidence by linear grating (1=oblique 0=collinear) 
B=1;    % expanding by fresnel lens (1=yes 0=no) 
L=7;   % topological charge 
Lambda=30;  % linear grating fringe spacing (in pixels) 
aspect=1/2;  % fraction of fringe period that is white/clear 
xsize=1024;  % resolution (xres x yres number data pts calculated) 
ysize=768;  % 

% define the X and Y ranges (defined to skip zero) 
xvec = linspace(-xsize/2, xsize/2, xsize);  % list of x values 
yvec = linspace(-ysize/2, ysize/2, ysize);  % list of y values 

% define the meshes - matrices linear in one dimension 
[xmesh, ymesh] = meshgrid(xvec, yvec); 

% calculate the individual phase components 
vortexPh = atan2(ymesh,xmesh);  % the vortex phase 
linPh = -2*pi*ymesh;   % a phase of linear grating 
radialPh = (xmesh.^2+ymesh.^2);  % a phase of defocus 

% combine the phases with appropriate scales (phases are additive) 
% the 'pi' at the end causes inversion of the pattern 
Ph = L*vortexPh + A*linPh/Lambda + B*radialPh/R^2; 

% transmittance function (the real part of exp(I*Ph)) 
T = cos(Ph); 
% the binary version 
binT = T > cos(pi*aspect); 

% plot the pattern 
% imagesc(binT) 
imagesc(T) 
colormap(gray)