## Double polar coordinates

The coordinate system for the auditory system.

We can hear the location of a sound source based on differences between two ears (binaural differences), and on the spectral filtering properties of the ear (spectral cues). Binaural differences are only useful for determining the horizontal location, and spectral cues mainly for the vertical and front-back locations. These two dimensions can best be described by the double-polar coordinate system.

I can’t seem to find any description of this coordinate system on Google or Wikipedia, which leads me to believe that it is not in wide practice. However, it is commonly used in scientific studies on sound localization: e.g. Knudsen and Konishi (1979) defined and used it in a study on the sound localization abilities of the barn owl.

# Double Polar Coordinate System

The double-polar coordinate system is a three-dimensional coordinate system in which each point on a sphere is determined by two angles (azimuth α and elevation ε) and a distance (depth r). The radial coordinate r denotes the point’s distance from a central point or the origin, and equals the radius of a sphere centered on the origin. The angular coordinates are defined (Fig. 1) relative to two reference planes through the origin: a horizontal and vertical plane, that are orthogonal to each other. Azimuth is the angular deviation from the vertical plane (and describes the horizontal deviation), while elevation is the angular deviation from the horizontal plane (and describes the vertical deviation).

On a unit sphere, the coordinate system looks composed of two orthogonal sets of circles, but appears as a rectangular, perpendicular grid of lines when viewed from the origin to $$\alpha , \epsilon$$=[0, 0] deg (Fig. 2).

The main advantage of such a double-polar coordinate system, is that the azimuth and elevation coordinates are independent of each other (unlike the spherical, Earth coordinate system). Another property of this coordinate system is that the absolute values of the azimuth and elevation coordinates of a single orientation in the frontal plane added cannot be larger than 90 deg: $$| \alpha | + | \epsilon | \leqq 90$$, so that all azimuth and elevation coordinates graphically all lie within a diamond shape (Fig. 3).

# Other coordinate systems

Double-polar coordinates can be converted to Cartesian coordinates by the trigonometric sine function and the Pythagorean Theorem:
$$x=r\sin(\alpha)$$
$$y=r\sin(\epsilon)$$
$$z=\sqrt{r^2-x^2-y^2}$$

Other coordinate systems include the spherical and polar coordinate system (Fig. 4).

The polar coordinate system is useful in describing rotations of the eye, while the spherical system is best known as a geographical coordinate system. The PANDA toolbox contains several m-files that can transform double-polar coordinates to these coordinates, and vice versa, such as:


pa_azel2cart; % transform double-polar to cartesian
pa_azel2pol; % transform double-polar to spherical
pa_rphi2azel; % polar to double-polar


# Why double-polar?

The cues for sound localization are interaural differences and spectral pinna cues, and these are oriented in (almost) the same fashion as the azimuth and elevation components of the double-polar coordinate system. This entails that if you perturb the interaural cues for horizontal sound localization (e.g. by plugging) or the spectral cues for vertical sound localization (e.g. by molding the pinna), only the azimuth or the elevation components will be independently affected. Using another coordinate system, both directions would appear to be affected.

This problem also holds for non-perturbed sound localization. Let's compare for example the double-polar with the spherical coordinate system. The spherical vertical coordinate $$\phi$$ is equivalent to the double-polar coordinate $$\epsilon$$. The horizontal spherical coordinate $$\theta$$ depends on both $$\alpha$$ and $$\epsilon$$: $$\theta = asin{\frac{sin{\alpha}}{cos{\epsilon}}}$$. Because localization is more imprecise in the elevation $$\epsilon$$ direction, this will be reflected also in the horizontal $$\theta$$ direction, thus seemingly introducing errors in sound localization (an example sound localization study is shown in Fig. 5).

This problem is not present when stimuli are presented and responses are made in the cardinal planes. Figure 6 shows how the double-polar $$\alpha$$ and the spherical $$\theta$$ coordinates deviate from each other, depending on the elevation $$\epsilon$$.

# Matlab Code

And finally, for those interested in creating their own figures, here is some not so well-documented Matlab code.

For figure 2:


function pa_doublepolar
%% Initialization
close all
clear all
clc

AZ  = -90:15:90;
AZ  = AZ';
n  = 100; % number of elevations at AZ(ii)
EL  = NaN(length(AZ),n); % create lots of potential elevations
for ii = 1:length(AZ)
mx = abs(round(abs(AZ(ii))-90));
if mx>0
y = linspace(-mx,mx,n);
else
y = zeros(1,n);
end
EL(ii,:) = y;
end
AZ    = repmat(AZ,1,n);
R    = 1;

%% Conversion
X       = R*sind(AZ); % Left-right
Y       = R*sind(EL); % Up-Down
signZ   = sign(cosd(AZ).*cosd(EL));
absZ    = abs(sqrt(R.^2-X.^2-Y.^2));
Z       = signZ .* absZ;

%% From the back to the front
% x,y,z scheme is rotated in Matlab, so we have to correct for this
X    = [X;fliplr(X)];
Y    = [Y;fliplr(Y)];
Z    = [Z;fliplr(-Z)];
AZ    = [AZ;fliplr(AZ)];
X    = X';
Y    = Y';
Z    = Z';
AZ    = AZ';

%% Graphics
for jj = [1 3 4]
subplot(2,2,jj);
plot3([-1 1],[0 0],[0 0],'b:','LineWidth',2);
hold on
plot3([-1 1],[0 0],[0 0],'bo','LineWidth',2,'MarkerFaceColor','b');
% let's give a color-gradient to azimuth
uAZ = unique(AZ(:)); % unique azimuths
nAZ = numel(uAZ); % number of unique azimuths
col = gray(nAZ); % a grayscale colormap
for ii = 1:nAZ
sel = AZ==uAZ(ii);
plot3(X(sel),Y(sel),Z(sel),'k-','LineWidth',2,'Color',col(ii,:));
end
axis([-1 1 -1 1 -1 1]);
axis square;
axis off;
end

%% Double-polar coordinates, now iso-elevation contours
EL  = -90:15:90;
EL  = EL';
AZ  = NaN(length(EL),n);
for ii = 1:length(EL)
mx = abs(round(abs(EL(ii))-90));
if mx>0
y = linspace(-mx,mx,n);
else
y = zeros(1,n);
end
AZ(ii,:) = y;
end
EL    = repmat(EL,1,n);
R    = 1;

%% Conversion
X       = R*sind(AZ); % Left-right
Y       = R*sind(EL); % Up-Down
signZ   = sign(cosd(AZ).*cosd(EL));
absZ    = abs(sqrt(R.^2-X.^2-Y.^2));
Z       = signZ .* absZ;

%% From the back to the front
X    = [X;fliplr(X)];
Y    = [Y;fliplr(Y)];
Z    = [Z;fliplr(-Z)];
EL    = [EL;fliplr(EL)];
% x,y,z scheme is rotated in Matlab, so we have to correct for this
X    = X';
Y    = Y';
Z    = Z';
EL    = EL';
[X,Y,Z] = pa_pitch(X(:),Y(:),Z(:),90); % just some rotation, see PANDA toolbox

%% Graphics
for jj = [2 3 4]
subplot(2,2,jj);
plot3([0 0],[0 0],[-1 1],'r:','LineWidth',2);
hold on
plot3([0 0],[0 0],[-1 1],'ro','LineWidth',2,'MarkerFaceColor','r');
% let's give a color-gradient to azimuth
uEL = unique(EL(:)); % unique azimuths
nEL = numel(uEL); % number of unique azimuths
col = gray(nEL); % a grayscale colormap
col = flipud(col);
for ii = 1:nEL
sel = EL==uEL(ii);
plot3(X(sel),Y(sel),Z(sel),'k-','LineWidth',2,'Color',col(ii,:));
end
axis([-1 1 -1 1 -1 1]);
axis square;
axis off;
end

subplot(2,2,4)
view(0,0);
axis([-1.01 1.01 -1.01 1.01 -1.01 1.01]);

%% Print
print('-depsc','-painter',mfilename); % I try to avoid bitmaps,
% so I can easily modify the figures later on in for example Illustrator
% For web-and Office-purposes, I later save images as PNG-files

function [X,Y,Z] = pa_pitch(X,Y,Z,Angle)
% [X,Y,Z] = PA_YAW(X,Y,Z,ROLL);
%
% Rotate YAW (deg) about the X-axis
%

% 2007 Marc van Wanrooij

%% Initialization
% Compensate for Coordinate Definition
[mx,nx] = size(X);
[my,ny] = size(Y);
[mz,nz] = size(Z);

% Compensate for Coordinate Definition in Matlab
x = X(:);
y = -Z(:);
z = Y(:);

%% Define Roll Rotation Matrix for 3D
R = [1 0          0           0;...
0  cos(Angle) -sin(Angle) 0;...
0  sin(Angle) cos(Angle)  0;...
0  0          0           1];

%% Define Data Matrix
M = [x y z ones(size(x))]';

%% Matrix Multiply
M = R*M;

%% Revert
M = M';
X = M(:,1);
Y = M(:,3);
Z = -M(:,2);

% Reshape
X = reshape(X,mx,nx);
Y = reshape(Y,my,ny);
Z = reshape(Z,mz,nz);


And here follows the code for figure 3.


%% Initialization
close all
clear all
clc

AZ  = -90:15:90;
AZ  = AZ';
n  = 100; % number of elevations at AZ(ii)
EL  = NaN(length(AZ),n); % create lots of potential elevations
for ii = 1:length(AZ)
mx = abs(round(abs(AZ(ii))-90));
if mx>0
y = linspace(-mx,mx,n);
else
y = zeros(1,n);
end
EL(ii,:) = y;
end
AZ    = repmat(AZ,1,n);
plot(AZ,EL,'.');

%% Graphics
for jj = [1 3]
subplot(1,3,jj);
% let's give a color-gradient to azimuth
uAZ = unique(AZ(:)); % unique azimuths
nAZ = numel(uAZ); % number of unique azimuths
col = cool(nAZ); % a cool colormap
for ii = 1:nAZ
sel = AZ==uAZ(ii);
plot(AZ(sel),EL(sel),'k-','LineWidth',2,'Color',col(ii,:));
hold on
end
end

%% Double-polar coordinates, now iso-elevation contours
EL  = -90:15:90;
EL  = EL';
AZ  = NaN(length(EL),n);
for ii = 1:length(EL)
mx = abs(round(abs(EL(ii))-90));
if mx>0
y = linspace(-mx,mx,n);
else
y = zeros(1,n);
end
AZ(ii,:) = y;
end
EL    = repmat(EL,1,n);

%% Graphics
for jj = [2 3]
subplot(1,3,jj);
hold on
uEL = unique(EL(:));
nEL = numel(uEL);
col = hot(nEL); % a hot colormap
col = flipud(col);
for ii = 1:nEL
sel = EL==uEL(ii);
plot(AZ(sel),EL(sel),'k-','LineWidth',2,'Color',col(ii,:));
end
end

for ii = 1:3
subplot(1,3,ii)
h = plot([-90 -0],[0 90],'k-',...
[0 90],[90 0],'k-',...
[0 90],[-90 0],'k-',...
[-90 0],[0 -90],'k-');
set(h,'LineWidth',2,'Color','k');
axis([-90 90 -90 90]);
axis square;
box off
set(gca,'XTick',-60:30:60,'YTick',-60:30:60);
xlabel('Azimuth (deg)');
ylabel('Elevation (deg)');
end

%% Print
print('-depsc','-painter',mfilename); % I try to avoid bitmaps,
% so I can easily modify the figures later on in for example Illustrator
% For web-and Office-purposes, I later save images as PNG-files


And this is the code to draw the sphere in figure 4.


function plotsph
TH    = 0:15:180;
PH    = -90:90;

[TH,PH] = meshgrid(TH,PH);
R    = ones(size(TH));
[X,Y,Z] = sph2cart(TH,PH,R);
for jj = 2
subplot(1,3,jj);
hold on
[m,n] = size(TH); % number of unique azimuths
col = flipud(gray(n)); % a grayscale colormap
for ii = 1:n
plot3(X(:,ii),Y(:,ii),Z(:,ii),'k-','LineWidth',2,'Color',col(ii,:));
end
axis([-1 1 -1 1 -1 1]);
axis square;
axis off;
end

TH    = -180:180;
PH    = -90:15:90;
[TH,PH] = meshgrid(TH,PH);
R    = ones(size(TH));
[X,Y,Z] = sph2cart(TH,PH,R);
for jj = 2
subplot(1,3,jj);
hold on
[m,n] = size(PH); % number of unique azimuths
col = flipud(gray(m)); % a grayscale colormap
for ii = 1:m
plot3(X(ii,:),Y(ii,:),Z(ii,:),'k-','LineWidth',2,'Color',col(ii,:));
end
axis square;
axis off;
end
subplot(1,3,2)
view(180,0);


And here is the code the draw the polar isocontour lines in figure 4.



function plotpol
TH    = -180:15:180;
PH    = -0:90;

[TH,PH] = meshgrid(TH,PH);
R    = ones(size(TH));
[X,Y,Z] = sph2cart(TH,PH,R);
for jj = 3
subplot(1,3,jj);
hold on
[m,n] = size(TH); % number of unique azimuths
col = flipud(gray(n)); % a grayscale colormap
for ii = 1:n
plot3(X(:,ii),Y(:,ii),Z(:,ii),'k-','LineWidth',2,'Color',col(ii,:));
end
axis([-1 1 -1 1 -1 1]);
axis square;
axis off;
end

TH    = -180:180;
PH    = -0:15:90;
[TH,PH] = meshgrid(TH,PH);
R    = ones(size(TH));
[X,Y,Z] = sph2cart(TH,PH,R);
for jj = 3
subplot(1,3,jj);
hold on
[m,n] = size(PH); % number of unique azimuths
col = flipud(gray(m)); % a grayscale colormap
for ii = 1:m
plot3(X(ii,:),Y(ii,:),Z(ii,:),'k-','LineWidth',2,'Color',col(ii,:));
end
axis([-1 1 -1 1 -1 1]);
axis square;
axis off;
end
subplot(1,3,3)
view(180,90);