Hi Lmac,
It isn't possible to have a coordinate system which has orthogonal axis and is shap preserving.
I came across a similar problem when trying to visualise the distribution of surface orientation for triangulation. My approach was to just normalise by the area instead. (little bit of a hack but had a good effect )
function [xs,ys,zs,as] = sphereHeatmap(az,el,w,n)
%% Function returns points for a heatmap mapped onto a unit sphere surface
% Inputs: az (n by 1) - azimuth (radians), el (n by 1) - elevation (radians)
% w (n by 1) - weighted value of points, n - number of bins is 2*n (integer)
% Outputs are of size (n+1) by (n+1)
% outputs can be used as surf(xs,ys,zs,as):
n_n = (-n:1:n);
%Round points to nearest integer n
azInt = round(az*2*n);
elInt = round(el*2*n);
%Generate spherical meshgrid
sintheta = sin(n_n*pi/n); sintheta(1) = 0; sintheta(2*n+1) = 0;
cosphi = cos(n_n'*pi/(2*n)); cosphi(1) = 0; cosphi(2*n+1) = 0;
xs = cosphi*cos(n_n*pi/n);
ys = cosphi*sintheta;
zs = sin(n_n'*pi/(2*n))*ones(1,2*n+1);
as = (cosphi*ones(1,2*n+1)); %Area of local small surface
indi = n_n==az;
indj = n_n==el;
for i = 1:2*n+1
for j = 1:2*n+1
if as(i,j)~=0
area_temp = sum(w((indi(:,i).*indj(:,j))==1),'all');
if as(i,j)~=0 % Ignore faces with zero surface area
as(i,j) = 100*area_temp./(sum(w)*(as(i,j))); %Normalise
end
end
end
end
end
