Auch wenn es schon mehr als drei Jahre her ist melde ich mich als Verfasser des ersten Beitrages auch noch einmal zu Wort.

Nach vielen Schwierigkeiten habe ich es irgendwann geschafft mein Projekt "Koordinatendatenbank auf Excelbasis" zu einem Ende zu bringen und nach kürzester Zeit war es aus dem Alltag unseres Vermessungsbüros nicht mehr wegzudenken, brauchte man doch nun nicht mehr das Katasteramt für Auskünfte anrufen etc.
Alles in allem ist es jedoch sehr speziell an die Ordner & Serverstruktur des Büros angepasst, so das es keinen Sinn macht dies hier zu veröffentlichen.
Jedoch scheint es nach wie vor Menschen zu geben, die das selbe Problem wie ich haben: Umrechnung und Transformation von Daten. Diesen Kern meines Programms auf Basis von Visual Basic for Applications (VBA) möchte ich gerne hier noch einmal veröffentlichen.

Vorraussetzung ist eine Koordinate mit einem Rechtswert = Variabel LKR und einem Hochwert = Variabel LKH

Für die Umrechnung einer UTM Koordinate in eine geografische Koordinate (in Dezimalform) gilt folgender Ansatz:

Public Function Transformation_UTMWGS84_nach_GeographischWGS84()
 
Const Pi = 3.14159265358979
m_zone = 32
 
'  WGS 84 Datum
'  Große Halbachse a
   a = 6378137#
'  Abplattung f
   F = 0.00335281068
 
'  Polkrümmungshalbmesser c
   C = a / (1 - F)
 
'  Quadrat der zweiten numerischen Exzentrizität
   ex2 = (2 * F - F * F) / ((1 - F) * (1 - F))
   ex4 = ex2 * ex2
   ex6 = ex4 * ex2
   ex8 = ex4 * ex4
 
'  Koeffizienten zur Berechnung der geographischen Breite aus gegebener
'  Meridianbogenlänge
   e0 = C * (Pi / 180) * (1 - 3 * ex2 / 4 + 45 * ex4 / 64 - 175 * ex6 / 256 + 11025 * ex8 / 16384)
   f2 = (180 / Pi) * (3 * ex2 / 8 - 3 * ex4 / 16 + 213 * ex6 / 2048 - 255 * ex8 / 4096)
   f4 = (180 / Pi) * (21 * ex4 / 256 - 21 * ex6 / 256 + 533 * ex8 / 8192)
   f6 = (180 / Pi) * (151 * ex6 / 6144 - 453 * ex8 / 12288)
 
 
'  Geographische Breite bf zur Meridianbogenlänge gf = LKH
   sigma = (LKH / 0.9996) / e0
   sigmr = sigma * Pi / 180
   bf = sigma + f2 * Sin(2 * sigmr) + f4 * Sin(4 * sigmr) + f6 * Sin(6 * sigmr)
 
'  Breite bf in Radianten
   br = bf * Pi / 180
   tan1 = Tan(br)
   tan2 = tan1 * tan1
   tan4 = tan2 * tan2
 
   cos1 = Cos(br)
   cos2 = cos1 * cos1
 
   etasq = ex2 * cos2
 
'  Querkrümmungshalbmesser nd
   nd = C / Sqr(1 + etasq)
   nd2 = nd * nd
   nd4 = nd2 * nd2
   nd6 = nd4 * nd2
   nd3 = nd2 * nd
   nd5 = nd4 * nd
 
'  Längendifferenz dl zum Bezugsmeridian lh
   lh = (Val(Mid(m_zone, 1, 2)) - 30) * 6 - 3
   dx = (LKR - 500000) / 0.9996
   dx2 = dx * dx
   dx4 = dx2 * dx2
   dx3 = dx2 * dx
   dx5 = dx4 * dx
   dx6 = dx3 * dx3
 
   b2 = -tan1 * (1 + etasq) / (2 * nd2)
   b4 = tan1 * (5 + 3 * tan2 + 6 * etasq * (1 - tan2)) / (24 * nd4)
   b6 = -tan1 * (61 + 90 * tan2 + 45 * tan4) / (720 * nd6)
 
   l1 = 1 / (nd * cos1)
   l3 = -(1 + 2 * tan2 + etasq) / (6 * nd3 * cos1)
   l5 = (5 + 28 * tan2 + 24 * tan4) / (120 * nd5 * cos1)
 
'  Geographischer Breite bb und Länge l als Funktion von Ostwert x und Nordwert y
   BB = bf + (180 / Pi) * (b2 * dx2 + b4 * dx4 + b6 * dx6)
   LL = lh + (180 / Pi) * (l1 * dx + l3 * dx3 + l5 * dx5)
 
End Function

Für die Umrechnung einer Gauß Krüger Koordinate in eine geografische Koordinate (in Dezimalform) brauche ich zwei Ansatz:

1. Ansatz:

Public Function Transformation_GaußKrüger_nach_GeographischBessel()
 
Const Pi = 3.14159265358979
 
'  Potsdam Datum
'  Große Halbachse a
   a = 6377397.155
'  Abplattung f
   F = 0.00334277321
 
'  Polkrümmungshalbmesser c
   C = a / (1 - F)
 
'  Quadrat der zweiten numerischen Exzentrizität
   ex2 = (2 * F - F * F) / ((1 - F) * (1 - F))
   ex4 = ex2 * ex2
   ex6 = ex4 * ex2
   ex8 = ex4 * ex4
 
'  Koeffizienten zur Berechnung der geographischen Breite aus gegebener Meridianbogenlänge
   e0 = C * (Pi / 180) * (1 - 3 * ex2 / 4 + 45 * ex4 / 64 - 175 * ex6 / 256 + 11025 * ex8 / 16384)
   f2 = (180 / Pi) * (3 * ex2 / 8 - 3 * ex4 / 16 + 213 * ex6 / 2048 - 255 * ex8 / 4096)
   f4 = (180 / Pi) * (21 * ex4 / 256 - 21 * ex6 / 256 + 533 * ex8 / 8192)
   f6 = (180 / Pi) * (151 * ex6 / 6144 - 453 * ex8 / 12288)
 
 
'  Geographische Breite bf zur Meridianbogenlänge gf = LKH
   sigma = LKH / e0
   sigmr = sigma * Pi / 180
   bf = sigma + f2 * Sin(2 * sigmr) + f4 * Sin(4 * sigmr) + f6 * Sin(6 * sigmr)
 
'  Breite bf in Radianten
   br = bf * Pi / 180
   tan1 = Tan(br)
   tan2 = tan1 * tan1
   tan4 = tan2 * tan2
 
   cos1 = Cos(br)
   cos2 = cos1 * cos1
 
   etasq = ex2 * cos2
 
'  Querkrümmungshalbmesser nd
   nd = C / Sqr(1 + etasq)
   nd2 = nd * nd
   nd4 = nd2 * nd2
   nd6 = nd4 * nd2
   nd3 = nd2 * nd
   nd5 = nd4 * nd
 
'  Längendifferenz dl zum Bezugsmeridian lh
   kz = Val(Mid(LKR, 1, 1))
   lh = kz * 3
   dy = LKR - (kz * 1000000# + 500000)
   dy2 = dy * dy
   dy4 = dy2 * dy2
   dy3 = dy2 * dy
   dy5 = dy4 * dy
   dy6 = dy3 * dy3
 
   b2 = -tan1 * (1 + etasq) / (2 * nd2)
   b4 = tan1 * (5 + 3 * tan2 + 6 * etasq * (1 - tan2)) / (24 * nd4)
   b6 = -tan1 * (61 + 90 * tan2 + 45 * tan4) / (720 * nd6)
 
   l1 = 1 / (nd * cos1)
   l3 = -(1 + 2 * tan2 + etasq) / (6 * nd3 * cos1)
   l5 = (5 + 28 * tan2 + 24 * tan4) / (120 * nd5 * cos1)
 
'  Geographischer Breite bp und Länge lp als Funktion von Rechts- und Hochwert
   BB = bf + (180 / Pi) * (b2 * dy2 + b4 * dy4 + b6 * dy6)
   LL = lh + (180 / Pi) * (l1 * dy + l3 * dy3 + l5 * dy5)
 
Call Transformation_GeographischBessel_nach_GeographischWGS84
 
End Function

2. Ansatz:

Public Function Transformation_GeographischBessel_nach_GeographischWGS84()
 
Const Pi = 3.14159265358979
 
'  Quellsystem Potsdam Datum
'  Große Halbachse a
   a = 6378137# - 739.845
'  Abplattung fq
   fq = 0.00335281066 - 0.00001003748
 
'  Zielsystem WGS84 Datum
'  Abplattung f
   F = 0.00335281066
 
'  Parameter für datum shift
   dx = 631     'Landesvermessung: 631 ; Ursprünglich 587
   dy = 23      'Landesvermessung: 23  ; Ursprünglich 16
   dz = 451     'Landesvermessung: 451 ; Ursprünglich 393

'  Quadrat der ersten numerischen Exzentrizität in Quell- und Zielsystem
   e2q = (2 * fq - fq * fq)
   e2 = (2 * F - F * F)
 
'  Breite und Länge in Radianten
   b1 = BB * (Pi / 180)
   l1 = LL * (Pi / 180)
 
'  Querkrümmungshalbmesser nd
   nd = a / Sqr(1 - e2q * Sin(b1) * Sin(b1))
 
'  Kartesische Koordinaten des Quellsystems Potsdam
   xp = nd * Cos(b1) * Cos(l1)
   yp = nd * Cos(b1) * Sin(l1)
   zp = (1 - e2q) * nd * Sin(b1)
 
'  Kartesische Koordinaten des Zielsystems (datum shift) WGS84
   X = xp + dx
   Y = yp + dy
   Z = zp + dz
 
'  Berechnung von Breite und Länge im Zielsystem
   rb = Sqr(X * X + Y * Y)
   BB = (180 / Pi) * Atn((Z / rb) / (1 - e2))
 
   If (X > 0) Then LL = (180 / Pi) * Atn(Y / X)
   If (X < 0 And Y > 0) Then LL = (180 / Pi) * Atn(Y / X) + 180
   If (X < 0 And Y < 0) Then LL = (180 / Pi) * Atn(Y / X) - 180
 
End Function

In allen Fällen kommt eine geografische Koordinate mit Längen- und Breitengrad heraus, deren Variablen LL und BB heißen.

Und nun wünsche ich allen viel Spaß beim tüfteln:
Schönen Gruß, Michael