Private Sub TRFtoGrid(Degrees As LatLon, a As Double, b As Double, N0 As Int, E0 As Int, F0 As Double, phi0 As Double, lambda0 As Double) As GridRef
' TRF Lat/Lon To Grid East/North by Transverse Mercator projection.
' <Phi> is latitude - negative if south of equator (all angles in radians)
' <Lambda> is longitude - negative if west of prime meridian
' <a> And <b> are major And minor axes of ellipsoid
' <N0> And <E0> are easting And northing of grid origin
' <F0> is scale factor on central meridian of grid
' <phi0> And <lambda0> are latitude And longitude of grid origin
' <East> And <North> are easting And northing of converted grid reference
' Tested And confirmed 11th Oct 2018
Dim Phi, Lambda As Double
Dim n, n2, n3, e2, nu, rho, eta2 As Double ' Geometric terms
Dim M, I, II, III, IIIa, IV, V, VI As Double ' Intermediate terms
Dim sin1, sin2, cos1, cos2, cos3, cos5 As Double
Dim tan1, tan2, tan4 As Double ' Trig's of <phi> and their powers
Dim phidiff, phi2diff, phi3diff As Double ' (Phi - phi0), (Phi + phi0) and
Dim phisum, phi2sum, phi3sum As Double ' their multiples
Dim L, L2, L3, L4, L5, L6 As Double ' (Lambda - lambda0) And its powers
Dim x As Double ' Temporary intermediate value term
Dim result As GridRef
Phi = Degrees.Latitude * RADIAN ' Convert to radian measure
Lambda = Degrees.Longitude * RADIAN
n = (a - b) / (a + b)
n2 = n * n
n3 = n2 * n
e2 = (a * a - b * b) / (a * a) ' Eccentricity squared value
' Calculate trigonometric factors
sin1 = Sin(Phi) : cos1 = Cos(Phi)
sin2 = sin1 * sin1
cos2 = cos1 * cos1 : cos3 = cos2 * cos1 : cos5 = cos3 * cos2
tan1 = sin1 / cos1 : tan2 = tan1 * tan1 : tan4 = tan2 * tan2
x = (1 - e2 * sin2)
nu = a * F0 / Sqrt(x)
rho = a * F0 * ( 1 - e2) / Sqrt(x * x * x)
eta2 = nu / rho - 1.0
phidiff = (Phi - phi0) : phi2diff = phidiff + phidiff
phi3diff = phi2diff + phidiff
phisum = (Phi + phi0) : phi2sum = phisum + phisum
phi3sum = phi2sum + phisum
M = b*F0 * ( (1 + n + 1.25 * (n2 + n3)) * phidiff _
- (3 * (n + n2) + 2.625 * n3) * Sin(phidiff) * Cos (phisum) _
+ 1.875 * (n2 + n3) * Sin(phi2diff) * Cos(phi2sum) _
- 35/24 * n3 * Sin(phi3diff) * Cos(phi3sum))
I = M + N0
x = nu * sin1
II = x * cos1 * 0.5
III = x * cos3 * (5 - tan2 + 9*eta2) / 24
IIIa = x * cos5 * (61 - 58*tan2 + tan4) / 720
IV = nu * cos1
V = nu * cos3 * (nu / rho - tan2) / 6
VI = nu * cos5 * (5 - 18*tan2 + tan4 + 14*eta2 - 58*tan2*eta2) / 120
L = (Lambda - lambda0)
L2 = L * L : L3 = L2 * L : L4 = L3 * L : L5 = L4 * L : L6 = L5 * L
result.Northing = I + II * L2 + III * L4 + IIIa * L6
result.Easting = E0 + IV * L + V * L3 + VI * L5
Return result
End SubDon't worry about getting my name wrong; it happens quite a lot and nobody would object to being called "Brain".
SELECT X(utmPt) AS X, Y(utmPt) AS Y
FROM (SELECT Transform(MakePoint([LON], [LAT], [SRID]), [outputSRID]) AS utmPt) AS Pt
you need a spatialite database in your project to make the calculations and dim it like
and the init function is
local_srid is the epsg code of yor system.
Hi George
[/
Sub PlotPoints(A As Int)
Dim AveY, AveX As Double
Dim cp As CameraPosition
Dim Co As CircleOptions
File.MakeDir(File.DirRootExternal & "/CEASER/","Data/Projections")
Co.Initialize
'Dim BitMap1 As Bitmap = LoadBitmapResize(File.DirRootExternal & "/CEASER/DATA/", "Point.png",8dip,8dip,True)
Try
gmap.Clear
If CGlobals.CoordCode=1 Then
'Site Coords
Dim rs As ResultSet = CGlobals.sql1.ExecQuery("SELECT Name, East, North, Elevation, Description FROM SCoords")
else if CGlobals.CoordCode=2 Then
'Global Coords
Dim rs As ResultSet = CGlobals.sql1.ExecQuery("SELECT Name, East, North, Elevation, Description FROM GCoords")
Else
'Job Coords
Dim rs As ResultSet = CGlobals.sql1.ExecQuery("SELECT Name, East, North, Elevation, Description FROM Coords")
End If
i=0: AveY=0: AveX=0
Do While rs.NextRow
Dim row(5) As Object
row(0) = rs.GetString("Name")
row(1) = NumberFormat2(rs.GetDouble("East"),1,3,3,False)
row(2) = NumberFormat2(rs.GetDouble("North"),1,3,3,False)
'row(3) = NumberFormat2(rs.GetDouble("Elevation"),1,3,3,False)
Geo.get_point_from_local(row(1),row(2))
Co.StrokeWidth(2)
' co.Center2(Engine.Lat,Engine.Lon).Radius(Radius.SelectedItem).FillColor(Colors.White).StrokeColor(Colors.DarkGray)
' gmapextra.AddCircle(gmap,co)
'gmap.AddMarker(Engine.Lat,Engine.Lon, "New Marker")
'gmap.AddMarker3(Engine.lat,Engine.Lon,row(0),Bitmap1)
i=i+1
AveY=AveY+row(1)
AveX=AveX+row(2)
Loop
rs.Close
If (AveY<>0 Or AveX<>0) And A=0 Then
'Zoom in to center of points
AveY=AveY/i
AveX=AveX/i
Convert(row(1),row(2),34)
' Engine.Hart94_YX(AveY,AveX,19)
' cp.Initialize(Engine.Lat,Engine.Lon, 10)
' gmap.MoveCamera(cp)
'ConvertLatLong
End If
Catch
Log(LastException)
End Try
End Sub
][/
Sub Class_Globals
Private Spa As Spatialite_Database
Dim PI As Double = 3.1415926535897932385
Dim RADIAN As Double = 2.0 * PI / 360
Dim Local_srid As String
'Dim Data_File As String
Dim spconstants As Spatialite_Constants
' Initialize(File.DirRootExternal & "/CEASER/Data/" ,Data_File,"READONLY")
' Spa.Open(File.DirInternal, "mm.sqlite", spconstants.SQLITE_OPEN_READONLY)
End Sub
Public Sub Initialize '(DPATH As String ,DBASE As String , OPEN_STYLE As String)
Dim DPath, DBase, Open_Style As String
Spa.Initialize
DPath=File.DirRootExternal & "/CEASER/Data"
DBase="Projections.sl3"
Open_Style="ReadOnly"
Spa.Open(DPath, DBase, spconstants.SQLITE_OPEN_READONLY)
Select Open_Style
Case "READWRITE"
Spa.Open(DPath, DBase, spconstants.SQLITE_OPEN_READWRITE)
Case "READONLY"
Spa.Open(DPath, DBase, spconstants.SQLITE_OPEN_READONLY)
End Select
End Sub
Sub get_point_from_local( x As Double , y As Double) As String
Dim sqStm As Spatialite_TableResult
Dim result As String
Dim sql1 As String
sql1 = "SELECT ST_AsText(ST_Transform(ST_GeomFromText('POINT("
sql1 = sql1 & x & " " & y & ")'," & CGlobals.EPSG & "),4326)) As trans_geom; "
Try
sqStm = Spa.GetTable(sql1)
result = sqStm.Rows(0,0)
Catch
Log(LastException.Message)
MsgboxAsync(LastException.Message & CRLF & "wrong EPSG code" ,"ATTENTION")
Local_srid="2048"
Return "0000000000000000000"
End Try
If result=Null Or result="" Then
Return "0000000000000000000"
End If
Return result
End Sub
'''transform lat lon to local
Sub get_point_to_Local( lat As Double , lon As Double) As String
Dim sqStm As Spatialite_TableResult
Dim result As String
Dim sql1 As String
sql1 = "SELECT ST_AsText(ST_Transform(ST_GeomFromText('POINT("
sql1 = sql1 & lon & " " & lat & ")',4326),"
sql1 = sql1 & CGlobals.EPSG & ")) As trans_geom;"
sqStm = Spa.GetTable(sql1)
result = sqStm.Rows(0,0)
Return result
End Sub
]
The EPSG code is correct, as I load it from my Database.
It keeps returning 0, when I want to convert from a Local System to Lats\Longs.
