Dear sig-geo listeRs,

I've been using an old version of the program geocalc (The geographic
calculator) to transform my field GPS data from the original WGS84 lat-long
to the local National_grid_system (Kandawala, in Sri Lanka)

I recently decided to update myself and try to use spTransform from the
library rgdal. It took me some time to figure out how to define Kandawala's
CRS. Once I that I (apparently) got it, I found that the resulting data show
a discrepancy with the results from geocalc. This difference is a constant
of 216 and 43 m for the x and y axis, respectively. I guess there might
something in my configuration.

I will appreciate if anybody can help me to identify the origin of this
difference or suggest a better way to proceed. Below is my script, with a
few points as example, plus some documentation.

Thank you very much,

Ahimsa


#=== start script

#=== data ====
# GPS points recorded in the field (Sri Lanka) using lat-long WGS84
ID <- c(1:4)
X <- c(81.19672, 81.19682, 81.19683, 81.19697)
Y <- c(6.20115, 6.20172, 6.20223, 6.20237)
wgs84.LL <- data.frame(ID,X,Y)

# the same data after transforming it using geocalc
X2 <- c(246811.2398, 246822.2554, 246823.3166, 246838.7971)
Y2 <- c(111592.7887, 111655.8285, 111712.2253, 111727.7191)
kand <- data.frame(X2,Y2)
# (below are the details of geocalc version and proj data)


#=== transformation process ===
# I create spatial data,
library(rgdal)
SP_wgs84.LL <- SpatialPoints(cbind(X,Y), proj4string=CRS("+proj=longlat
+datum=WGS84"))

# I define Sri Lankan National grid CRS,
my.kandawala <- "+proj=tmerc +ellps=evrst30 +lat_0=7.000480277777778 +lon_0=
80.771711111111130 +k=0.9999238418 +x_0=200000.000 +y_0=200000.00"

# and then I change the projection:
k.points <- spTransform(SP_wgs84.LL, CRS(test))
k.points

#=== comparing the points from geocalc and spTransform:
kand - as.data.frame(k.points)
# there is a difference of ~216m and ~43m on the x and y axis, respectively

#============== end script

Documentation

About Geocalc:
I'm using Geocalc, the geographic calculator version 3.05 (1992-1994). The
original points in degrees (dd.ddddd) and configured as "system=Geodetic",
"Datum=186 WGS 1984" are converted into a user defined Kandawala Datum. The
details of this datum are stored in a DAT file. These are the contents in
such file:

1000 Zone_1 0 4 0 4 85 0 200000.000000000000000 200000.000000000000000
80.771711111111130 7.000480277777778 0.999923840000000


# Kandawala in ArcGIS
ArcGIS contains a prj file named kandawala. These are its contents:

GEOGCS["GCS_Kandawala",DATUM["D_Kandawala",SPHEROID["Everest_1830",
6377276.345,300.8017]],PRIMEM["Greenwich",0],UNIT["Degree",
0.017453292519943295]]

# I found more details about Kandawala on the web:

http://www.mail-archive.com/mapinfo-l@csn.net/msg01850.html

# kandawala data

#     FUGRO SURVEY - GEO Version 2.38.01

#      GEODESY 1
#      Datum                 :  Kandawala
#      Spheroid              :  Everest 1830C
#      Semi Major Axis       :  6377276.345 m
#      Inverse Flattening    :  300.801700000
#      Projection            :  Transverse Mercator (UTM)
#      Latitude  Origin      :    7ų 00' 01.7290" N
#      Longitude Origin      :   80ų 46' 18.1600" E
#      False Easting         :       200000.000 m
#      False Northing        :       200000.000 m
#      Central Scale Factor  :  0.9999238418


# PARAMETERS FOR CONVERSION FROM WGS 84
#                               Kandawala      WGS 84       Geod 1 to Geod 2
#      dX                    :   98.300   m      0.000   m    -98.300   m
#      dY                    : -787.500   m      0.000   m    787.500   m
#      dZ                    :  -85.000   m      0.000   m     85.000   m
#      rX                    :    0.0000  "      0.0000  "      0.0000  "
#      rY                    :    0.0000  "      0.0000  "      0.0000  "
#      rZ                    :    0.0000  "      0.0000  "      0.0000  "
#      dS                    :    1.00000 ppm    0.00000 ppm   -1.00000 ppm

#      TRANSFORMATIONS
#      Station   Name        :  PW-5
#                Datum       :  Kandawala            :  WGS 84
#                Projection  :  Transverse Mercator  :  Transverse Mercator
#                Latitude    :    6ų 57' 09.3865" N  :    6ų 57' 10.5800" N
#                Longitude   :   79ų 50' 37.0926" E  :   79ų 50' 44.7687" E
#                Height      :            8.097 m    :          -93.704 m
#                Easting     :        97459.221 m    :       372492.852 m
#                Northing    :       194807.060 m    :       768701.916 m
#                PSF         :     1.0000539617      :     0.9998012206
#                Convergence : -  0ų 06' 44.47"      :  -  0ų 08' 23.08"





-- 
ahimsa campos-arceiz
www.camposarceiz.com

	[[alternative HTML version deleted]]


