Close

Star Observation Code

A project log for Daylight Geolocator Remix

This is a remix of https://hackaday.io/project/28550-light-level-geolocator

agpcooperagp.cooper 03/22/2018 at 01:294 Comments

In the Beginning

In the beginning I was going to use "star observation" code but I was attracked to the direct sunrise/sunset "latitude" code. Why not, it is what is presented as the "solution" on all the Internet sites I visited.

The star obsevation code is more general and is based on a graphical iterative method. I wrote some code for my programmable calculator 35 years ago for it.

The problem with the latitude code is that it assumes the Declination and Right Ascension are constant. They are not, so it is not very accurate. It can be fixed but it is not a very tidy solution (just average the noon and midnight estimates).

Here is my version of the latitude and longitude codes:

double latitude(double SunRise,double SunSet,double Alt,double Dec) {
  double DL,A,B,C,D,U,V,Lat;

  DL=SunSet-SunRise;
  A=tan(radians(Dec));
  B=cos(radians(7.5*DL));
  C=sin(radians(Alt))/cos(radians(Dec));
  D=sqrt(A*A+B*B);
  if (fabs(C)<=D) {
    U=degrees(asin(B/D));
    V=degrees(asin(C/D));
    Lat=V;
    if (DL>0) Lat=V-U;
    if (DL<0) Lat=V+U;
    if (Dec<0) Lat=-Lat;
  } else {
    Lat=0.0;
  }
  return Lat;
}


And:

    JD=J2000(Year,Month,Day);
    SunData(JD,(Sunrise+Sunset)/2.0-LTC,&RA,&Dec);
    Lat=latitude(Sunrise,Sunset,Alt,Dec);
    GST=GMST(JD,(Sunrise+Sunset)/2.0-LTC);
    Lon=15.0*(RA-GST);
    if (Sunrise>Sunset) Lon=fmod(Lon-180,180);


If you use this code with exact sunrise and sunset data you get:

And:

The star observation code is more general:

void StarObs(float RA1,float Dec1,float Alt1,float UT1,float RA2,float Dec2,float Alt2,float UT2,float *Lat,float *Lon)
{
  float P1,P2,U1,U2,V1,V2,X0,Y0,LHA1,LHA2,EstLat,EstLon;
  EstLat=radians(*Lat);
  EstLon=radians(*Lon);
  RA1=radians(RA1*15.0);
  Dec1=radians(Dec1);
  Alt1=radians(Alt1);
  LHA1=radians(UT1*15.0);
  RA2=radians(RA2*15.0);
  Dec2=radians(Dec2);
  Alt2=radians(Alt2);
  LHA2=radians(UT2*15.0);
  do {
    P1=asin(sin(Dec1)*sin(EstLat)+cos(Dec1)*cos(EstLat)*cos(EstLon+LHA1-RA1));
    P2=asin(sin(Dec2)*sin(EstLat)+cos(Dec2)*cos(EstLat)*cos(EstLon+LHA2-RA2));
    U1=(P1-Alt1)*cos(Dec1)*sin(EstLon+LHA1-RA1)/cos(P1);
    U2=(P2-Alt2)*cos(Dec2)*sin(EstLon+LHA2-RA2)/cos(P2);
    V1=(Alt1-P1)*(sin(Dec1)-sin(EstLat)*sin(P1))/cos(EstLat)/cos(P1);
    V2=(Alt2-P2)*(sin(Dec2)-sin(EstLat)*sin(P2))/cos(EstLat)/cos(P2);
    if (V1*U2!=U1*V2) {
      X0=((V1*V1+U1*U1)*U2-(V2*V2+U2*U2)*U1)/(V1*U2-U1*V2);
      Y0=((V2*V2+U2*U2)*V1-(V1*V1+U1*U1)*V2)/(V1*U2-U1*V2);
    } else {
      X0=0.0;
      Y0=0.0;
    }
    *Lat=*Lat+X0;
    *Lon=*Lon+Y0;

  } while ((fabs(X0)>0.000001)||(fabs(Y0)>0.000001));
  *Lat=degrees(EstLat);
  *Lon=degrees(EstLon);

Now there are two questions about the code. Is it robust and how good does the initial estimate need to be?

In my head (based on the algorithm) the initial estimates need to be within 90 degrees (i.e latitude). The reason is the code assumes the earth is flat so 90 degrees is basically infinity. 

Is it robust? Tests say no! It did not work for altitudes less than about -4 degrees. It seems fine for the range of interest (i.e. -3 to 0 degrees):

The above plot show the star observation results:

Note the variance (Blue line) is zero for -0.833 degrees.

This plot is not quite what it seems. The data is exact for -0.833 degrees and the variance reflects the use of an incorrect altitude. The algorithm appears stable for altitudes "estimates" from -3 to 0 degrees.

From the plot the latitude is very sensitive to the altitude estimate.

The plot also indicates that for good qualitity data it is possible to reverse engineer the observations altitude.

The Longitude estimate hardly varies.

For uncorrected (for light levels) experimental data gave the following results:

Searching for the minimum variance gave:

Which is off 3 degrees in latitude (almost exact in longitude) but improvements are possible if lighting conditions are taken into account.

AlanX

Discussions

agp.cooper wrote 03/26/2018 at 03:37 point

Hi Ted,

I will have a look. At the moment my task is to get the sunrise/sunset intensity correction factor working. Had a "bomb" last week collecting data (coding error somewhere). At the moment I can estimate Lat, Lon and Alt if I have multiple readings and the altitude is consistent. That is the algorithm determines the altitude of the observations as well!

So yeah, I can use a camera to look for the sun (or moon) and determine latitude and longitude when it crosses and arbitary horizontal line. For the stars I need to know the actual altitude as the RA and Dec do not change enough to be measurable.

The maths behind star postions and star observations are not processing intensive. I am using a variation of Newton's method and it takes about 5 interations. Each iteration has five points (calculations) for the estimating f' and f". Finding/matching a star goup with an arbitary scale, direction and rotation would be hard (OpenCV?). Star observations will require some known reference point(s) such as an altitude or an azimuth, however.

However, using a "levelled" camera to watch the sun (and the moon) seems like a really good follow up project. I have the maths to solve the problem.

AlanX

  Are you sure? yes | no

agp.cooper wrote 03/25/2018 at 08:16 point

Ted,

Long long time ago before the days of GPS I used an almanac, a stopwatch and a radio to workout where I was (when lost in the outback of WA looking for gold). Never found much gold (a few colours) but I didn't get lost for long (only overnight).

We sighted (measured the time and altitude) of two stars. You can do the same with a good clinometer or sextant.

Did you know some of the brightest stars  (e.g. Sirus) are visible to the naked eye in daylight (if you know where to look)!

Anyway, I never thought I would ever use the code again.
AlanX

  Are you sure? yes | no

Ted Yapo wrote 03/25/2018 at 17:28 point

Yes, I've seen a few stars, Venus, and Iridium flares during daylight.  It's a fun thing to try.

Oh, and the Moon :-)

There are systems that determine attitude of satellites or high-flying aircraft from sighting stars like this; it's a very similar problem.  Modern systems can recognize the stars automatically.  You need a camera and a bunch of processing power, but it's probably more precise than tracking sunrise/sunset.  There is a large body of literature on the subject.  This survey paper is a bit dated, but covers a lot of the basics:

http://www.mdpi.com/1999-4893/2/1/93/htm

  Are you sure? yes | no

Ted Yapo wrote 03/22/2018 at 01:50 point

Star observation?  Are you working on the lost-in-space problem, or manually sighting one star?

  Are you sure? yes | no