3

I have a solar panel pointing (it's normal vector) in some direction. I want to calculate the angle between that and the current position of the sun. I am using pyephem and I have this information in two pairs of azimuth and altitude. 

panel_az = ephem.degrees('180')
panel_alt = ephem.degrees('45')
sun_az = ephem.degrees('245')
sun_alt = ephem.degrees('22')

What is the easiest way to find the angle between the panel's normal vector and the vector pointing towards the sun?

Mads Skjern
  • 5,648
  • 6
  • 36
  • 40

3 Answers3

5

The library offers a separation() function that gives the angle between two spherical coordinates; look near the bottom of this section of the Quick Reference:

http://rhodesmill.org/pyephem/quick.html#other-functions

I think you will get the angle you are seeking if you run:

a = ephem.separation((panel_az, panel_alt), (sun_az, sun_alt))
print a

Good luck!

Brandon Rhodes
  • 83,755
  • 16
  • 106
  • 147
2

Convert both to vectors first:

z = sin(altitude)
hyp = cos(altitude)
y = hyp*cos(azimuth)
x = hyp*sin(azimuth)
vector = (x,y,z)

Then calculate the angle between the vectors (say a and b) using cross and dot products.

angle = atan2(norm(cross(a,b)), dot(a,b))

For cross use:

def cross(a, b):
    c = [a[1]*b[2] - a[2]*b[1],
         a[2]*b[0] - a[0]*b[2],
         a[0]*b[1] - a[1]*b[0]]
    return c

For dot use:

def dot(a, b):
    c = [ a[i] * b[i] for i in range(len(a)) ]
    return c

For norm use:

def norm(a):
    mag = sqrt(sum(a[i]*a[i] for i in range(len(a))))
    c = [ a[i]/mag  for i in range(len(a)) ]
    return c
Henk van der Laak
  • 305
  • 2
  • 7
  • I guess that works. I did smth similar before. However, I was hoping for a solution that would not require so much coding from my side. PyEphem is big, so I was hoping this could be calculated with some convenience functions. Also, SciPy is huge, so is it really necessary to write my own definition of cross, dot and norm? – Mads Skjern Sep 08 '13 at 20:02
0
clear,clc,clf;
for n=[15 46 74 105 135 166 196 227 258 288 319 349];
delta=23.45*(sind(360*(284+n)/365));
phi=31.2;
omega=acosd(-tand(phi).*tand(delta));
omega1=-omega:0.1:omega;
alt=(sind(delta).*sind(phi))+(cosd(delta).*cosd(omega1).*cosd(phi));
alt1=asind(abs(alt))
azm=(cosd(delta).*sind(omega1))./cosd(alt);
azm1=asind(azm)
 plot(azm1,alt1,'b');hold on;
  end
grid on
xlabel('Solar Azimuth');ylabel('Solar altitude');
text(0,85,'noon')
hold off
Mohamed Elsayed
  • 1