function v1_mapping
% V1_MAPPING  Illustrate afferent mapping in V1.
%
%   V1_MAPPING() plots a 2-panel figure:
%     (1) Visual space (retinal error): fixation F and target T
%     (3) v1 neural space: complex-log map
%
%

%% Initialisation
% clear;
clc;
close all;



%% IsoR and isoPhi contour grids
figure(1);
clf
%% Iso-R-contours
R		= [1 2.3 5.4];
phi		= -180:180;
[R,phi] = meshgrid(R,phi);
phi		= phi/180*pi;
[x,y] = pol2cart(R,phi);
x		= R.*cos(phi);
y		= R.*sin(phi);
subplot(1,3,1)
plot(x,y,'k-','LineWidth',1.5,'Color',[.8 .8 .8]);
hold on
nicegraph;
axis equal;
set(gca,'XTick',[],'YTick',[])
title('visual space');

R		= [1 2.3 5.4];
nR = numel(R);
for ii = 1:nR
[x,y] = pol2cart(R(ii),45);
str = ['$r = ' num2str(R(ii)) '^\circ$'];
text(x,y,str,HorizontalAlignment="center",VerticalAlignment="bottom",FontSize=16,Interpreter="latex",Rotation=-45);
end

%% Iso-phi-contours
R		= linspace(0,5.4,100);
phi		= -180:45:180;
[R,phi] = meshgrid(R,phi);
phi		= phi/180*pi;
x		= R.*cos(phi);
y		= R.*sin(phi);

plot(x',y','k-','LineWidth',1.5,'Color',[.8 .8 .8]);
nicegraph;
xlim([-7 7]);
ylim([-7 7]);

%% Isocontours for neural space
R			= [1 2.3 5.4];
phi			= linspace(-90,90,1000);
[R,phi]		= meshgrid(R,phi);
[u,v]		= v1map(R,phi);

subplot(133)
plot(u,v,'k-','LineWidth',1.5,'Color',[.8 .8 .8]);
hold on

R			= linspace(0,5.4,1000);
phi			= -90:45:90;
[R,phi]		= meshgrid(R,phi);
[u,v]		= v1map(R,phi);

plot(u',v','k-','LineWidth',1.5,'Color',[.8 .8 .8]);
nicegraph;
xlabel('X (mm)','FontSize',21);
ylabel('Y (mm)','FontSize',21);
title('neural space (V1)');
xlim([-1 25]);
ylim([-20 20]);

R		= [1 2.3 5.4];
nR = numel(R);
for ii = 1:nR
[u,v] = v1map(R(ii),-45);
str = ['$r = ' num2str(R(ii)) '^\circ$'];
text(u,v,str,HorizontalAlignment="center",VerticalAlignment="bottom",FontSize=16,Interpreter="latex",Rotation=-90);
end
%% Save
fname = fullfile('..','..','images','v1_mapping_ori.svg');
savegraph(fname,'svg');

function [X,Y] = v1map(R,phi)
% V1MAP  Complex-logarithmic V1 motor map.

l	= 12; %mm
R0	= 1; % deg


X = l*log(1+R/R0);
Y = -l*R.*phi*pi./(180*(R0+R));

end


function [x,y] = pol2cart(R,phi)
% POL2CARTD  Polar-to-cartesian conversion with angles in degrees.

x		= R.*cosd(phi);
y		= R.*sind(phi);
end