[Proj] How to build the correct paramater list for pj_init()

[Roger Oberholtzer]

On Mon, 2007-10-29 at 12:57 -0600, Richard Zinn wrote:
> Hi, Sorry if this has been addressed before, I couldn't see any way to
> search the archive.
>  
> I want to build the correct paramater list for proj based on this
> data:
> PROJCS["NAD_1983_UTM_Zone_12N",GEOGCS["GCS_North_American_1983",
> DATUM["D_North_American_1983",
> SPHEROID["GRS_1980",6378137.0,298.257222101]],
> PRIMEM["Greenwich",0.0],
> UNIT["Degree",0.0174532925199433]],
>  
> PROJECTION["Transverse_Mercator"],
> PARAMETER["False_Easting",500000.0],
> PARAMETER["False_Northing",0.0],
> PARAMETER["Central_Meridian",-111.0],
> PARAMETER["Scale_Factor",0.9996],
> PARAMETER["Latitude_Of_Origin",0.0],
> UNIT["Meter",1.0]]
>  
> I've read the documentation, and since I'm just a newbie to
> cartography it is all still a little unclear and I need to learn by
> example.
>  
> Is there any code out there (any language) I could see that takes this
> kind of input and builds out the paramaters for prj_init()?
>  
> My guess is that it is something like this:
> static char *params_in[] = { 
>   "proj=tmerc", 
>   "lat_0=0", 
>   "lon_0=500000", 
>   "units=m" 
> };
> prj_init(4, params_in);
>  
> But that doesn't work :)  And I'm sure it is pretty far off.

Here is an example I use in the UK. I use a different routine, but you
should be able to do the same:        

Initialize thus: 

	projPJ fromProj, toProj;

	fromProj = pj_init_plus(
						"+proj=latlong "
						"+ellps=WGS84 "
						"+towgs84=0,0,0 "
						"+nodefs");

	toProj   = pj_init_plus(
						"+proj=tmerc "
						"+lat_0=49 "
						"+lon_0=-2 "
						"+k=0.9996012717 "
						"+x_0=400000 "
						"+y_0=-100000 "
						"+ellps=airy "
						"+nadgrids=gb.gsb "
						"+units=m "
						"+no_defs");

fromProf describes the source data. In my case, WGS84 from a Trimble
receiver. toProj describes what I want out. In my case, OSTN02 in the
UK. At least this is how I interpret things. This code works for me.

To do a conversion, LONG is the longitude and LAT is the latitude, both
in degrees. I do not do anything with the altitude (HEIGHT):

	#define DEGREE_TO_RADIAN (M_PI / 180.0)

	lng = LONG * DEGREE_TO_RADIAN;
	lat = LAT  * DEGREE_TO_RADIAN;
	alt = 0; // We do not give this a height

	pj_transform(fromProj, toProj, 1, 0, &lng, &lat, &alt);

	local_easting = lng;
	local_northing = lat;
	local_altitude = HEIGHT;

The tricky part is the parameters to use!


>  
> Thanks.
> _______________________________________________
> Proj mailing list
> Proj at lists.maptools.org
> http://lists.maptools.org/mailman/listinfo/proj
-- 
Roger Oberholtzer

OPQ Systems / Ramböll RST

Ramböll Sverige AB
Kapellgränd 7
P.O. Box 4205
SE-102 65 Stockholm, Sweden

Tel: Int +46 8-615 60 20
Mobl: Int +46 70-815 1696