Смещение пикселя в тайле

[EXVIK]

Всем привет!
Не ругайтесь, я начинающий, я пока трудно вникаю во все это...
Интересно вот что.
Для сфероида смещение пикселя получается вычислять верно.
Мне было сказано:
Точно также как ты вычисляешь путь к картинке только зум считать на 8 больше.
Получишь координаты пиксела в огромной картинке.
Дальше уже пересчитывай эти координаты пиксела в экранные координаты. Ну и касательно тайла тоже...
Вот код чутка:

Код: Выделить всё

NumTiles:=Power(2,Round(Rec.Zoom)+8-1);
BitmapSize:=NumTiles * 256;
BitmapOrigo:= BitmapSize / 2;
PixelsPerLonDegree := BitmapSize / 360;
PixelsPerLonRadian := BitmapSize / (2*PI);
X:=Floor(BitmapOrigo + Rec.Lon * PixelsPerLonDegree);
Z:=Sin(Rec.Lat*(pi/180));
Y:=Floor(BitmapOrigo - 0.5 * ln((1+z)/(1-z)) * PixelsPerLonRadian);
X:=Floor(X/256);            //X на большой картинке
Y:=Floor(Y/256);            //Y на большой картинке
X2:=Floor(X/256);           //В каком блоке по X
Y2:=Floor(Y/256);           //В каком блоке по Y
ShiftPixel.X:=X-X2*256;   //Смещение X
ShiftPixel.Y:=Y-Y2*256;   //Смещение Y


Но на эллипсоид это смещение не верно!
Подскажите, каким образом высчитывать смещение для эллипсоидной проекции?
Отмечу, что Y вычислял не так, а вот так:

Код: Выделить всё

Z:=Sin(Rec.Lat*(pi/180));
Y:=Floor(BitmapOrigo - PixelsPerLonRadian*(Arctanh(Z)-Exct*Arctanh(Exct*Z)));


Т. е. заменял 2 строчки, но результат не верный. Помогите разобраться...

[EXVIK]

Переформулирую свою проблему.
Зная точно Lat, Lon и Zoom, как мне найти точно пиксель, которому будет соответствовать эти данные? Для Сфероидальной проекции, вроде, вычисляет верно(вычисление выше), а для Эллипсоидальной проекции чутка промахивается для моего города это не значительно...2-3 пикселя при Zoom=16. Как правильно считать пиксель, в который указывают точные данный Lat, Lon, Zoom?

[snake32]

Присоединяюсь к вопрсу.
Какую формулу использовать для перевода LatLonToXY и обратно для тайлов яндекса?
Для гугла и яндекса формулы разные как я понял.
Формулу гугла нашёл:

Код: Выделить всё

procedure LonLatToXY( const Lon,Lat:Extended; zoom:integer; var x,y:Extended );
   var z:Extended;
begin
  z := sin (Lat*toRad);
  X := BitmapOrigo[zoom] + Lon * PixelsPerLonDegree[zoom];
  Y := BitmapOrigo[zoom] - 0.5 * ln((1+z)/(1-z)) * PixelsPerLonRadian[zoom];
end;

Вывел обратную:

procedure XYToLonLat( const x,y:Extended; zoom:integer; var Lon,Lat:Extended );
var n:Extended;
begin
Lon := (x - BitmapOrigo[zoom]) / PixelsPerLonDegree[zoom];
n := exp(-2*(y-BitmapOrigo[zoom])/PixelsPerLonRadian[zoom]);
Lat := arcsin( (n-1)/(n+1) )*toDeg;
end;

А для яндекса тайлов уже вторую неделю ищу. Подскажите пожалуйста, люди добрые!!!
Ссылку http://ru.wikipedia.org/wiki/%D0%9F%D1% ... 1%80%D0%B0
понять не могу. Если можно, кусочек кода с переводом.

[snake32]

Все формулы которые удалось загуглить возвращают неверный результат:

Код: Выделить всё

{
  double rLat = lat * Math.PI / 180;
  double a = 6378137.0;
  double k = 0.0818191908426;

  double zz = Math.Tan(Math.PI / 4 + rLat / 2)  / Math.Pow((Math.Tan(Math.PI / 4 + Math.Asin(k * Math.Sin(rLat)) / 2)), k);
  double y = (20037508.342789 - a * Math.Log(zz)) * 53.5865938 / Math.Pow( 2 ,23 - z);
  return ((long)y );
  или
      z = Tan(pii / 4 + rLat / 2)/(Tan(pii / 4 + ArcSin(k * Sin(rLat)) / 2)) ^ k
      GetPix.X = CLng((20037508.342789 +a * rLon) * 53.5865938 / 2 ^ (23 - Mstb))
      GetPix.Y = CLng((20037508.342789 -a * Log(z)) * 53.5865938 / 2 ^ (23 - Mstb))
}

function LatToYwgs84( const Lat:Extended; zoom:Cardinal ):Extended;
  var n,rLat:Extended;
begin
  rLat := Lat*toRad;
  n := power( tan( pid4+rLat*0.5 ) / tan( pid4 + arcsin( k*sin( rLat ) )*0.5), k );
  Result := ( 20037508.342789 - a*Ln(n) )*53.5865938/Power( 2 ,23-zoom );
end;

     {double a = 6378137.0;
            double c1 = 0.00335655146887969;
            double c2 = 0.00000657187271079536;
            double c3 = 0.00000001764564338702;
            double c4 = 0.00000000005328478445;


            double mercY = 20037508.342789 - (y * Math.Pow(2,23 - z)) / 53.5865938;
            double g = Math.PI / 2 - 2 * Math.Atan(1 / Math.Exp(mercY / a));
            double zz = g + c1 * Math.Sin(2 * g) + c2 * Math.Sin(4 * g) + c3 * Math.Sin(6 * g) + c4 * Math.Sin(8 * g);
            return (zz * 180 / Math.PI);}
function YtoLatWGS84( const y:Extended; zoom:Cardinal ):Extended;
  var mercY,g,zz:Extended;
begin
  mercY := 20037508.342789 - (y * Power(2,23 - zoom) ) / 53.5865938;
  g := pid2 - 2 * arctan(1 / exp(mercY / a) );
  zz := g + c1*sin(2*g) + c2*sin(4*g) + c3*sin(6*g) + c4*sin(8*g);
  Result := zz*toDeg;
end;

I need HELP!!!