Hess coordinate transformation (sphere projection) I'm running into an issue with a specific coordinate transformation. 
What I got so far
The projection is the following: A sphere with center (0, 0, 0) and Radius 1. On this sphere, there are arc segments, always with 5° separation:

Behind the sphere is a plane, that intersects the sphere at a single point and is perpendicular to the normal axis. Now we trace the segment lines from the center point (0, 0, 0) and project them onto the plane. If we do that, we get this pattern:

This is my coordinate system. Every 5 degrees is one step in the system (see annotations). I can draw the below pattern (the screenshot comes from my program with the points added with photoshop) like this: I draw the horizontal lines by calculating the (x, y) like this:

These get me the traces when following the red line (chose theta and then let phi run). The blues I get by rotating the whole thing, doing x = y, y = -x. 
Where I'm stuck
This transformation doesn't get me the points where the horizontal and the vertical lines meet, because the crossection will non-linearely 'distort' away. In the end, I want to draw in this coordinate system, also with coordinates like (3.4, -2.3) and stuff like that. Also, I need to be able to project back from screen coordinates to the original coordinates, so I need a good transformation that gives me the exact position.
What is the correct way to do this? If there are polynominal approximations, that would work too, it doesn't need to be absolutely mathematically correct, just have this shape. Also, this view is the only one we need, we don't need different lines, it'a always -45° to 45°.
I'm probably going about this all wrong, but this distortion effect keeps throwing me off...
(The chart that I need to implement is called 'Hess Chart" and is used in an opthalmologic sub-speciality called 'orthoptics' (mostly diagnosis and treatment of cross-eyed children)).
Thanks in advance!
EDIT: It looks like if we would flatten a paraboloid, we would get a similar pattern. 

 A: Starting from $x=\tan\phi$, $y=\cot\theta\sec\phi$, where $\phi$ is free and $45°\le\theta\le135°$ because $\theta=90°$ corresponds to $y=0$ we can recover $\theta$ because
$$\frac{y^2}{x^2+1}=\frac{\cot^2\theta\sec^2\phi}{\tan^2\phi+1}=\cot^2\theta$$
And so
$$\cot\theta=\frac{y}{\sqrt{x^2+1}}$$
Because $y\ge0$ for $\theta$ in the first quadrant.  
Similarly the vertical stripes are given by $x=\cot\alpha\sec\beta$ with $\beta$ free and $45°\le\alpha\le135°$. We can solve this for
$$\cot\alpha=\frac x{\sqrt{y^2+1}}$$
This gets us from $(x,y)$ on the screen to $(\alpha,\theta)$ on the sphere.  
To get to the crossing points, $y^2\tan^2\theta=\sec^2\phi=\tan^2\phi+1=x^2+1$. Seeing what $x$ and $y$ are in terms of $\alpha$ and $\beta$,
$$\tan^2\beta\tan^2\theta=\cot^2\alpha\sec^2\beta+1=\cot^2\alpha\left(\tan^2\beta+1\right)+1=\cot^2\alpha\tan^2\beta+\csc^2\alpha$$
So
$$\begin{align}\csc^2\alpha&=(\tan^2\theta-\cot^2\alpha)\tan^2\beta=(\sec^2\theta-\csc^2\alpha)\tan^2\beta\\
\tan^2\beta&=\frac{\csc^2\alpha}{\sec^2\theta-\csc^2\alpha}=\frac{\cos^2\theta}{\sin^2\alpha-\cos^2\theta}\\
\tan\beta&=\frac{\cos\theta}{\sqrt{\sin^2\alpha-\cos^2\theta}}\end{align}$$
Because both sides have the same sign as $y$. Similarly,
$$\tan\phi=\frac{\cos\alpha}{\sqrt{\sin^2\theta-\cos^2\alpha}}$$
So knowing $\theta$ and $\alpha$ we can find $\phi$, hence $x$ and $y$ on the screen.  
As an example we plot the grid as projected onto the screen in Matlab. Also we show how to transform a point from $(\alpha,\theta)$ coordinates to screen coordinates.
% Hess.m

clear all;
close all;

phi = linspace(-45,45,300);
for theta = 50:5:130,
    x = tand(phi);
    y = 1/tand(theta)./cosd(phi);
    if mod(theta,15) == 0,
        plot(x,y,'r-');
    else
        plot(x,y,'r:');
    end
    hold on;
end

beta = linspace(-45,45,300);
for alpha = 50:5:130,
    x = 1/tand(alpha)./cosd(beta);
    y = tand(beta);
    if mod(alpha,15) == 0,
        plot(x,y,'b-');
    else
        plot(x,y,'b:');
    end
    hold on;
end

axis([-1 1 -1 1]);
axis square;

xs = -5;
ys = -6;
alpha = 90-xs*5;
theta = 90-ys*5;
phi = atand(cosd(alpha)/sqrt(sind(theta)^2-cosd(alpha)^2));
x = tand(phi);
y = 1/(tand(theta)* cosd(phi));
plot(x,y,'kx');
text(x+0.05,y,['(' num2str(xs) ',' num2str(ys) ')']);
title('Projection Onto Screen');
xlabel('x');
ylabel('y');


