onekiloparsec / SwiftAA

The most comprehensive collection of accurate astronomical algorithms in (C++, Objective-C and) Swift.
http://www.onekiloparsec.dev/
MIT License
172 stars 31 forks source link

Fix localMidnight #54

Closed andris-zalitis closed 7 years ago

andris-zalitis commented 7 years ago

In some cases localMidnight(longitude:) returns wrong value.

let jd = JulianDay(year: 2016, month: 12, day: 20, hour: 3, minute: 5, second: 3.5)

let longitude4 = 90.0.degrees
print(jd.localMidnight(longitude: longitude4).date)
// Prints 2016-12-20 06:00:00 +0000

If we're on longitude W 90, and current time is December 20th, 03:05 UTC, then the date is not actually December 20th for us yet. December 20th will come when time is December 20th, 06:00 UTC.

onekiloparsec commented 7 years ago

Hello Andris. Thanks for questioning this. It is really a tricky question, and we may need to actually change the method signature instead.

This method is intended to not refer to any date object, as it is a complex concept. The idea is to have a simple mean to have a geographic/longitudinal shift from midnight UTC. In that sense, it is not the midnight in the usual human sense. Move left or right in longitude by some kilometers, and the value will change, whatever the timezone.

In that sense, the method provides this convenient shift. Mut maybe we should find some better naming?

andris-zalitis commented 7 years ago

I see. I have say I had a suspicion that it could be that this is the intended way. So it's only about naming then.

Perhaps we could have it like this:


public func localMidnight(for longitude: Degree) -> JulianDay 
public func midnight(for timeZone: TimeZone) -> JulianDay
andris-zalitis commented 7 years ago

Or perhaps the new one could be

public func midnightOfCurrentDate(for timeZone: TimeZone) -> JulianDay
andris-zalitis commented 7 years ago

I was just thinking. Is localMidnight really doing what you want it to do?

Look at this use case:

func twilights(forSunAltitude sunAltitude: Degree, coordinates: GeographicCoordinates) -> (rise: JulianDay?, transit: JulianDay?, set: JulianDay?, error: CelestialBodyTransitError?) {

    var error: CelestialBodyTransitError? = nil

    let jd_midnight = self.julianDay.localMidnight(longitude: coordinates.longitude)
    let sun_midnight = Sun(julianDay: jd_midnight)
    if sun_midnight.makeHorizontalCoordinates(with: coordinates).altitude > sunAltitude {
        error = .alwaysAboveAltitude
    }
    ...

Given Sun instance on a particular JulianDay, you convert the JulianDay to local midnight to check whether Sun actually goes below requested altitude (the rationale being that at midnight it might be approximately at the lowest point).

OK, let's think about it. Say we:

The date that we created corresponds to 2017-01-01 07:00:00 +11h (Sydney time)

So what will localMidnight do?

So when calculating twilights for date that corresponds to local time 2017-01-01 07:00:00 we will check whether sun was below altitude at midnight for previous day i.e. 2016-12-31 01:55:12. I think this will give us incorrect results because we should have checked Sun altitude at 2017-01-01 01:55:12 (2016-12-31 13:55:12 +0000) instead.

andris-zalitis commented 7 years ago

There is another, perhaps simpler perspective.

So the final result may differ from the original date by almost -1.5 JD. That will not be the nearest local midnight.

onekiloparsec commented 7 years ago

Hm, that's a thorough check of this problematic midnight. I am not the original author of this function. If I read the associated documentation, it reads

Returns the Julian Day corresponding to the Greenwhich midnight before the actual value

But the first -0.5 operation is there probably to ensure somehow to round correctly towards a date before the actual value. As see it today, this comes with an additional .rounded(.down). This makes no sense to have both.

I think that this method could simply be written return value.rounded(.down) and it would satisfy the definition, and fix the problems you describe here above, isn't it?

andris-zalitis commented 7 years ago

I think midnight is correct.

public var midnight: JulianDay { return JulianDay((value - 0.5).rounded(.down) + 0.5) }

Both 0.5 values are there because Julian Days are counted from noon time instead of midnight.

Since Julian Days are counted from noon, midnight in UTC will always have 0.5 at the end. That's why if we want to round down hours after the midnight, we remove 0.5, round and then add 0.5 back.

Having return value.rounded(.down) would return noon time in UTC before the actual value.

onekiloparsec commented 7 years ago

You're right. I am trying to play with the code right now, and I see I'm mistaken.

So what is wrong is to consider that, to compute midnight in somewhere else than Greenwich, one needs to compute it first for Greenwich and then shift it by longitude, no?

So we keep midnight as it is (maybe renaming it as midnightUTC or something else), and work on making localMidnight(for longitude|timeZone) correct.

onekiloparsec commented 7 years ago

While playing with the code, I think the handling of errors in the function func twilights(forSunAltitude sunAltitude: Degree, coordinates: GeographicCoordinates) -> (rise: JulianDay?, transit: JulianDay?, set: JulianDay?, error: CelestialBodyTransitError?) is wrong too.

I mean, I use computations outside the RiseTransitSetTimes function to infer values .alwaysAbove|BelowAltitude. That's risky. I'll work to better infer meaningful errors from basic error flags returned by RiseTransitSetTimes.

andris-zalitis commented 7 years ago

So what is wrong is to consider that, to compute midnight in somewhere else than Greenwich, one needs to compute it first for Greenwich and then shift it by longitude, no?

Yes.

So we keep midnight as it is (maybe renaming it as midnightUTC or something else), and work on making localMidnight(for longitude|timeZone) correct.

Well,localMidnight(for timeZone:) is already correct. At least as far as I know (time zones are tricky :))

public func midnight(for timeZone: TimeZone) -> JulianDay {
    let offsetFromGMT = JulianDay(Double(timeZone.secondsFromGMT(for: self.date)) / (60*60*24))
    return (self + offsetFromGMT).midnight - offsetFromGMT
}

midnight is calculated for offset value, and then the offset is removed.

onekiloparsec commented 7 years ago

I've juste committed a new implementation of localMidnigt(for longitude) method, covered by additional tests. I think, I got it right now. Can you have a look, and possibly close the issue? Thanks!

P.S. #55 becomes obsolete, isn't it?

andris-zalitis commented 7 years ago

I will take a look. But I wonder why didn't you use the one-liner I introduced in #55? :)

public func localMidnight(longitude: Degree) -> JulianDay {
    return (self - longitude.inHours.inJulianDays).midnight + longitude.inHours.inJulianDays
}

Does the new implementation returns different results? Also please note that we'd still need localMidnight for time zone because of #56

onekiloparsec commented 7 years ago

Interestingly, there is one single test that fails with your one-liner. I couldn't investigate more, but here is a screenshot. Interesting.

screen shot 2017-03-16 at 11 07 42
andris-zalitis commented 7 years ago

Ah, interesting. I will check this.

andris-zalitis commented 7 years ago

Oh, take a look at this:

public struct Degree: NumericType, CustomStringConvertible {
      ...
      public var inHours: Hour { return Hour(value / 15.0) }
      ...
}

let longitude8 = -360.0.degrees
print(longitude8.inHours)
// prints: -24h0m00.000s
// is this correct?
onekiloparsec commented 7 years ago

Well, yes this is technically correct, since there is no automatic "reduction" to a range like, say, [-180,180[ for instance. But I am unsure of where does this question lead us somewhere?

Note that I've worked quite a bit on the rise transit and set functions. This is because I wanted to reproduce the values produced by the C++ tests. Unfortunately, these tests only print values out and dont really test anything. When trying to reproduce the rise and set times of Mercury in March 1st, 2038 (as in the C++ tests), I fixed quite a few details. But I also broke some other ones. And if you look at the code today, some tests do not pass anymore.

andris-zalitis commented 7 years ago

But I am unsure of where does this question lead us somewhere?

I asked because I used inHours in my one-liner for localMidnight and thought that it provides reduced hours for a given longitude. Let's throw away my implementation and use yours 👍🏻

Note that I've worked quite a bit on the rise transit and set functions. This is because I wanted to reproduce the values produced by the C++ tests.

Looks like you don't use localMidnight anymore for RiseTransitSet. I'm OK with that, I will pre-calculate midnights for time zone myself. That actually means that we don't need midnight(for: TimeZone) anymore.

When trying to reproduce the rise and set times of Mercury in March 1st, 2038 (as in the C++ tests), I fixed quite a few details. But I also broke some other ones. And if you look at the code today, some tests do not pass anymore.

I see, can't comment much regarding those tests. I'm mostly checking Sun and Moon data against USNO.

onekiloparsec commented 7 years ago

A quick note. I have not checked more rise transit and set stuff for Sun and Moon, but I managed to get solar system objects times right. I'm not using localMidnight, as AA book is very precisely stating that coordinates must be computed for 0h Dynamical Time (which is basically midnight UT, not the local midnight). I'll thus close this issue.