Where an eclipse and crescent can be seen on a map?

[ebraminio]

Hi! Wondering if is possible to see where an eclipse and/or crescent will be visible?

See for example https://www.moonsighting.com/visibility.html and https://svs.gsfc.nasa.gov/4906

Initially my first question was how to draw day and night map which the question led me to http://thenewkinetic.blogspot.com/2011/05/how-to-draw-daynight-world-map_7282.html

Thus I have the following very initial result on my foss app :)

[cosinekitty]

This is an interesting problem. Reading your link to "The New Kinetic", I would use a similar idea. You need some way to calculate the latitude and longitude for every pixel on your map. It looks like you are using a Mercator projection. But this would work for any projection, so long as you have a way of calculating the geographic coordinates for a given pixel.

First of all, there are solar eclipses and lunar eclipses. Both will require different approaches. A lunar eclipse and a crescent moon observation will both be very similar, because they are asking "can you see the Moon from this place?" The solar eclipse is asking, "is the Moon blocking sunlight from reaching this place?"

I will have to think more about the solar eclipse case. The closest thing I have right now is SearchLocalSolarEclipse, which will tell you the beginning, middle, and end of a solar eclipse at a given location on the Earth. The problem is you would have to call it for every single pixel on your map, which would be extremely slow. Maybe I can think about this to find a better answer.

Lunar eclipses and crescent moons are more feasible with the existing Astronomy Engine code. You will want to find the moment in time you are interested in for observing the Moon. For example, perhaps you call SearchLunarEclipse to find a lunar eclipse. This will tell you when the lunar eclipse starts, peaks, and ends. Let's assume we want to make a map for a single moment in time, such as the peak of a lunar eclipse.

Calculating the Moon's position is expensive (the formulas are very complicated). So you will want to calculate its position just once for your entire map. The best way would be to call GeoMoon to obtain the Moon's geocentric position vector at whatever moment in time you are interested in.

Then, for each pixel on your map, calculate the geographic coordinates (latitude, longitude). Use those coordinates to create an Observer object that you can pass into a call to ObserverVector to calculate a geocentric vector of that point on the Earth's surface. Subtract the observer vector from the moon vector to obtain a vector from the observer toward the moon.

Then call Rotation_EQJ_HOR to calculate a rotation matrix that converts equatorial J2000 vectors into local horizon vectors, for that same observer location. (This part could use some optimization to avoid redundantly calculating the Earth's axis and rotation for every single pixel on your map. I will have to come back and think about that some more.)

Pass your vector from the observer toward the moon, and the rotation matrix, into a call to RotateVector. The result will be a horizon-based vector for the given location on the Earth. The x component points north, the y component points west, and the z component points up. If z is a positive number, it means the object is above the observer's horizon, and is thus visible.

This might be enough to get you started. Again, this is very interesting! I may consider adding some new code to Astronomy Engine to make map calculations like this more efficient, because I can see how the approach I outlined above will waste time by repeating a lot of calculations for every pixel on the map.

[ebraminio]

Thank you so much! Just to note the map is from here https://commons.wikimedia.org/wiki/File:World_map_blank_without_borders.svg and is equirectangular projection AFAIK, and using https://jakearchibald.github.io/svgomg/ and Zopfli for compression it only takes 36,498 bytes of gzip, I wished it had the zero longitude on its center which apparently isn't and I don't know where exactly it starts from, am playing with QGIS right now to see if I can generate another map maybe also, just said as you may want to do the same here or not :)

[cosinekitty]

@ebraminio Which programming language would you want to use to generate your maps? I was assuming JavaScript, but I thought I would ask, so that I can prototype something in that language first. I want to focus on eliminating the redundant Earth axis calculations for every pixel.

[ebraminio]

Oh, thanks a lot for asking, feel free to use any of the supported languages you feel will be easier for it :)

[ebraminio]

Guess drawing day night also would be great for example I applied one simple change on my codes persian-calendar/persian-calendar@f3667ff which made day night calculations much faster

[cosinekitty]

That is interesting also, because it turns out to be the same problem as the lunar eclipse case, only for the Sun instead of the Moon. In a way, it will be easier because, unlike for the Moon, there is no need for parallax correction for the Sun, because it is so much farther away. I will keep this in mind as I design the new feature.

[ebraminio]

Just to update, this is the result of my work about day night map till now (does't share any code with astronomy engine but guess still is interesting)

Screen.Recording.2022-01-20.at.11.50.13.PM.mov

It uses a equirectangular projection map I crafted https://commons.wikimedia.org/wiki/File:Equirectangular_projection_world_map_without_borders.svg using QGIS which only takes 50kb of gzip.

am letting you know so possibly you can reuse some of the ideas of it 😊

Update:

Now features moon also!

[ebraminio]

Sorry for repeated spam again :) day night map probably can be done faster also using https://astronomy.stackexchange.com/a/43818 but as I've only have spotted this haven't yet studied it yet in detail, I will try to implement it in my code base to see what experiences I can share, maybe something similar can be done for other maps in question.

[cosinekitty]

No, it's not spam at all; it's helpful! That is a good point, that at the resolution of a world map, it may not make any visible difference to treat the Earth as a perfect sphere by ignoring its oblateness. I'm still planning on adding a feature for enabling more efficient repeated calculations that could be used for many different applications, not just drawing maps, so I don't want to assume a spherical Earth in what I'm doing. However, it's great to point out other reference materials here for anybody else who is reading this, so they can make the best choice for their own project.

[cosinekitty]

I just pushed this commit: 90a9839. It improves performance for making maps like we have been discussing. I realized the main source of wasted calculation was the internal calculations for the Earth's rotation at a given moment in time, also known as sidereal time. The optimization is to cache the sidereal time and recycle it as much as possible. I did this for all 4 supported programming languages.

In the C version, for example, a simulated map-making program runs about 2.3 times faster with the sidereal time caching optimization.

I was very happy that the changes were simpler than I thought, and they do not require any changes in the exported interface. The changes are hidden entirely inside Astronomy Engine without breaking any outside code.

I'm going to convert this issue into a discussion, so we can keep it around for future reference, but without considering it an open issue.

[ebraminio]

Excellent! Thanks :)

[cosinekitty]

Here is a C++ example function that generates a Mercator projection of the Earth in a PNG image like below. Yellow shows intensity of sunlight, and blue shows intensity of moonlight. This program finishes in less than 0.8 seconds on my system. The image size is 1800 x 800 pixels, for a total of 1,440,000 pixels. Hopefully you can use a lower resolution and make this much faster.

[cosinekitty]

OK, now you have AstroVector.toObserver and Observer.toVector. They are inverse functions, and they correctly handle the fact that the Earth bulges at the equator.

[cosinekitty]

Also, you don't need to wait on RotateVector, because you already have RotationMatrix.rotate.

[ebraminio]

Both fixed now https://github.com/persian-calendar/persian-calendar/blob/34ffe54c701e0bc822092bacc9b48df535570a70/PersianCalendar/src/main/kotlin/com/byagowi/persiancalendar/ui/map/MapScreen.kt#L235

Let me know if anything else needs to be improved. Thank you so much! :)

[cosinekitty]

That is great! You could simplify verticalComponents a little bit:

private fun verticalComponent(rot: RotationMatrix, ovec: AstroVector, bvec: AstroVector): Double {
    val hor = rot.rotate(bvec - ovec)
    return hor.z / hor.magnitude()
}

Also, instead of searching for sunX, sunY, moonX, moonY inside the loops, I think you could calculate them once outside the loops. But I'm not sure which way is faster. It might not make any measurable difference. In my sample C program, I show how to find the geographic location directly beneath a celestial body, if you want to try it. You would use AstroVector.toObserver for this.

[ebraminio]

Applied https://github.com/persian-calendar/persian-calendar/blob/c8c0921a2272525417caf2593a5a2d8276508b49/PersianCalendar/src/main/kotlin/com/byagowi/persiancalendar/ui/map/MapScreen.kt#L210 Thanks! :)

But I'm not sure which way is faster. It might not make any measurable difference.

Felt slower AFAICT, guess because the current approach isn't making less than one degree accuracy but this tries to generate an accurate one.

[cosinekitty]

@ebraminio : I created a function checklist spreadsheet to help track which functions have been ported to Kotlin, so I can make sure I don't miss anything. I will keep this updated as I work. You might find it a handy reference.

[ebraminio]

Fantastic, I will have an immediate for use those eclipses I already have https://github.com/persian-calendar/persian-calendar/blob/073d859aeaea2d668dde279ca35346ece4fda705/PersianCalendar/src/main/kotlin/com/byagowi/persiancalendar/ui/astronomy/Eclipse.kt which I can replace with that hopefully.

[cosinekitty]

OK, I will give the eclipse calculations priority.

[ebraminio]

I now use astronomy's eclipse logic over my own, btw you might find this interesting :)

Screen.Recording.2022-04-18.at.4.38.48.AM.mov

[ebraminio]

One thing that can be very helpful is this map,

brought from https://www.moonsighting.com/visibility.html

Do you know how those curves are drawn? I guess it is something based on both moon altitude and elongation/azimuth. The map in general very accurately predicts where new moon crescent can be seen which I guess is an interesting topic by its own.

[ebraminio]

Just as a report maybe, so far we were able to recreate one of the maps @hidp123 and I wanted accurately https://github.com/crescent-moon-visibility/crescent-moon-visibility/blob/main/visibility.cc and the fact astronomy-engine provides both Python and C port has been very good for us as @hidp123 is more active on figuring out what should be done in the astronomical algorithms side using Python APIs and I'm more active on its other sides in C/C++.

[ebraminio]

Also turned out magnetic declination and other magnetic related calculations are already provided in Android https://developer.android.com/reference/android/hardware/GeomagneticField and the following is drawn using one of the calculations it provides

and here is its source

https://cs.android.com/android/platform/superproject/+/master:frameworks/base/core/java/android/hardware/GeomagneticField.java

[cosinekitty]

Wow, I did not realize the Earth's magnetic field was that complicated!

[ebraminio]

True, probably something was wrong there, here is what I ended up with which apparently is more correct, hopefully you would like it also:

Screen.Recording.2022-09-12.at.2.22.49.AM.mov

[ebraminio]

Thanks to @hidp123 helps now I've reached to this initial work,

Screen.Recording.2022-09-16.at.6.48.33.PM.mov

This needs more work https://github.com/persian-calendar/persian-calendar/blob/1d273d21598d624e35bd2545b81ba8e2b7caddbc/PersianCalendar/src/main/kotlin/com/byagowi/persiancalendar/ui/map/MapScreen.kt#L427 specially I think it has many duplicated calculations but now I wish I can learn how to draw lunar and solar eclipses visibility maps also in future as this at least have reached to somewhere 😊