Why are conversion results from proj.4 inaccurate?

Question

I want to convert a Gauss-Kruger coordinate (EPSG::31467) into a WGS84 coordinate (EPSG::4326) using proj.4. The problem is that the output from proj.4 is very, very inaccurate.

Consider the Gauss-Kruger coordinate (3509502, 5399989). I converted it with three different online converters. All results agreed to within four decimal places: (48.7373, 9.1281). Moreover, I verified with Google Maps that (48.7373, 9.1281) is indeed the correct WGS84 coordinate (I know where it should be).

However, proj.4 converts the Gauss-Kruger coordinate into (48.7384, 9.1292), more than 140 meters away from the correct solution. Why?

$ echo 3509502 5399989 | invproj +proj=tmerc +lat_0=0 +lon_0=9 +k=1 +x_0=3500000 +y_0=0 +ellps=bessel +datum=potsdam +units=m +no_defs -f "%.4f"
9.1292  48.7384

The results are the same for version 4.8.0 and the newest version 4.9.3.

---

The solution (from @JGH):

$ echo 3509502 5399989 | cs2cs -f %.12f \
> +proj=tmerc +lat_0=0 +lon_0=9 +k=1 +x_0=3500000 +y_0=0 +datum=potsdam +no_defs \
> +to +proj=longlat +datum=WGS84 +no_defs
9.128151599064  48.737368376923 51.846968705766

Answer 1

From the fine documentation:

The proj program is limited to converting between geographic and projection coordinates within one datum.

and

Output geographic coordinates will be in DMS and precise to 0.001" with trailing, zero-valued minute-second fields deleted.

The task you are trying is changing datum (Postdam 83 to WGS 84)

Answer 2

The epsg file in my GDAL directory is using these definitions

# DHDN / 3-degree Gauss-Kruger zone 3
<31467> +proj=tmerc +lat_0=0 +lon_0=9 +k=1 +x_0=3500000 +y_0=0 +datum=potsdam +units=m +no_defs  <>

Gdaltransform utility http://www.gdal.org/gdaltransform.html gives results which should make your happy

gdaltransform -s_srs epsg:31467 -t_srs epsg:4326
3509502 5399989
9.12815159906357 48.737368376923 51.8469687057659

Because also gdaltransform is using proj4 I believe that problem in not in proj4 library itself but in the way you are using it. However, I can't say what is wrong or what is missing in your command if it is not the +ellps=bessel that is different in your proj string.