There's nothing wrong with the formula, but it uses radians throughout. So you should convert your "latitude = 50.0" to "latitude = 0.87266..." as the point where the formula passes 1.
I need to draw maps on a computer screen.
From Wikipedia, I found that for a given latitude, the "y" value for a Mercator projection is given by...
y = ln(tan(Latitude) + (1 / Cos(Latitude)))
... so I wrote a program to calculate "y" for latitudes from 0.0 to 89.9.
Now from a "sit back and look at it" position, the "y" value increases the further you move from the equator, so if I have a 1 degree square cell, when I plot it, with say 100 pixels for 1 degree of longitude, the scaling due to the projection should mean that I multiply that same 100 in longitude by the factor I get from the formula above, and that I should end up with a figure larger then 100, i.e. that "y" calculated above should be greater then 1 for all values of latitude.
It is not. I have to set latitude = 50.0 before the result of the above formula passes 1. I am multiplying the latitude in degrees by Pi/180.0 to convert the degrees to Radians. I am using 3.14159265358979323846 as Pi. All of my variables are either declared to be double precision floating points, or explicitly stated with a decoimal point to force float operation by the compiler.
Is the formula wrong? Am I doing something silly?
I am multiplying the latitude in degrees by Pi/180.0 to convert the degrees to Radians. I am using 3.14159265358979323846 as Pi.
When my loop variable is 50.0 it is converted to Radians by multiplying by Pi/180.0 giving 0.87266462599716477 which goes into the formula and comes out as 1.0106831886830212.
That is even worse. When I do that, the value returned by the formula is less then 1 for all values.
The Wikipedia article says...
I read that as this being a formula which scales Latitude and Longitude to x and y coordinates on a screen, (or whatever).The following equations determine the x and y coordinates of a point on a Mercator map from its latitude φ and longitude λ (with λ0 being the longitude in the center of map):
This is the inverse of the Gudermannian function:
This is the Gudermannian function:
The scale is proportional to the secant of the latitude φ, getting arbitrarily large near the poles, where φ = ±90°. Moreover, as seen from the formulas, the pole's y is plus or minus infinity.
What I see is that, assuming the centre of the map is zero, x = Longitude - zero, so I draw the point 00:00:00N 00:00:00E at the origin. The point 00:00:00N 01:00:00E at (1 - 0) * (how many pixels per degree). and so on. Latitude being Zero.
The point 01:00:00N 00:00:00E is still drawn at x = 0, (from above), but the y point is now...
rLat = 1.0 * (180.0 / Pi)
y = ln(tan(rLat) + (1 / cos(rLat))) * (how many pixels per degree)
... to generate the necessary exageration as I move away from the equator.
To acheive the exageration y must surely be >1 for all non-zero values, but it is not.
Okay, let's use that example of Latitude = 50°. If then and (just working to 2 decimal places). Then my formula for y gives . That is somewhat greater than 50, as one would wish.
In fact, the "exaggeration" in y at latitiude is equal to the derivative of , which is , and this is always greater than 1.
Arrgh!
When I do that with a calculator, it works as you say, however, the program I have here does not give the same answer. The 1.19 and 1.56 agree, I can see them in the debugger, but then the next step gives a different result.
I will look at it further.
Fixed!
Works now, you won't believe how stupid I was being, when you said convert the degrees to Radians, I knew I had done that, but also knew I had "the conversion factor", already calculated as I didn't want the overhead of calculating it every time.
Trouble is, the factor I had was Pi / 180.0 which is fine for where I was using it, but not where I added it later on your reco!
Suffice to say, I now have two factors r2d and d2r.
I hope I don't need another r2d factor otherwise I will have r2d2 and risk Imperial intervention.