GPS to UTM

I have to translate GPS WGS84 coord Lat\Lon on UTM-WGS84 coord,
I found on internet a procedure written in c++ but it doesn't work
properly
if I compare the output with the output of other translator software,
and I don't know why. Maybe because I'm in Italy and I have to setup
different
parameter value?

How can I solve the problem? Where can I find a documentation or a
tested procedure?
Any one can help me?






deg2rad := PI / 180.0;

//WGS-84. See list with file "LatLong- UTM conversion.cpp" for id
numbers
//ReferenceEllipsoid WGS84
a := 6378137.000; //WGS84
ellipsoid[ReferenceEllipsoid].EquatorialRadius;
        e2 := 0.00669438; //WGS-84 //
ellipsoid[ReferenceEllipsoid].eccentricitySquared;
        k0 := 0.9996; //riduzione di scala

//Make sure the longitude is between -180.00 .. 179.9
        LongTemp := (Long+180) - Trunc((Long+180)/360)*360-180; // -180.00 ..
179.9;

        LatRad := Lat * Deg2Rad;
        LongRad := LongTemp * Deg2Rad;

        ZoneNumber := Trunc((LongTemp + 180)/6) + 1; //Trunc = integer

(LongTemp < 12.0) ) then
                ZoneNumber := 32;

// Special zones for Svalbard
        begin
31
ZoneNumber := 33
ZoneNumber := 35
ZoneNumber := 37;
        end;
        LongOrigin := (ZoneNumber - 1)*6 - 180 + 3; //+3 puts origin in
middle of zone
        LongOriginRad := LongOrigin * deg2rad;

        //compute the UTM Zone from the latitude and longitude
        UTMZone := Format('%d%s', [ZoneNumber, UTMLetterDesignator(Lat)]);

        eps := (e2)/(1-e2);

        N := a/sqrt(1-e2*sin(LatRad)*sin(LatRad));
        T := tan(LatRad)*tan(LatRad);
        C := eps*cos(LatRad)*cos(LatRad);
        A := cos(LatRad)*(LongRad-LongOriginRad);

        M := a*((1 - e2/4        - 3*e2*e2/64 - 5*e2*e2*e2/256)*LatRad
                                - (3*e2/8        + 3*e2*e2/32        + 45*e2*e2*e2/1024)*sin(2*LatRad)
                                + (15*e2*e2/256 + 45*e2*e2*e2/1024)*sin(4*LatRad)
                                - (35*e2*e2*e2/3072)*sin(6*LatRad));

        UTMEasting :=
(k0*N*(A + (1-T+C)*A*A*A/6 + (5-18*T+T*T+72*C-58*eps)*A*A*A*A*A/
120))
+ 500000.0;

        UTMNorthing := k0*(M + N * tan(LatRad)*
(A*A/2+(5-T+9*C+4*C*C)*A*A*A*A/24 + (61-58*T+T*T
+600*C-330*eps)*A*A*A*A*A*A/720));

        if (Lat < 0) then
                UTMNorthing := UTMNorthing + 10000000.0; //10000000 meter offset for
southern hemisphere

Gianluca wrote:

Perhaps some of these can be used to check results:
http://edu-observatory.org/maps/utm.html
http://edu-observatory.org/maps/utm_grid.html

Gianluca wrote:

Gianluca,

I think I can help you. The code that was in your post looks like
Pascal. I have a source code file of many pascal routines (Delphi7) that
I wrote to do a whole lot of map related stuff. In there is my function
LL2UTM (latitude/longitude to UTM) taken from the book "Map Projections
- A Working Manual" by John P. Snyder. You are welcome to the source
code to use as you please. Please understand that this is NOT an
application but rather a set of routines you can build your own
application. There is also the function UTM2LL that goes the other way.
Included in the code are is the function VInverse which is Vincenty's
Inverse Geodetic Solution and VDirect which is his Direct Geodetic Solution.

I also have the same code written in Visual Basic 6 as well. Let me know
if you are interested. You can email me directly at:
windswaytoo@gmail.com

Please put [MAPSTUFF] in the subject line.

Hope this helps

Vic

Thank you everybody for the help, I finally found my mistake!!
When I translate c++ code I didn't pay attention that there was two
variable with the same name a and A (wonderfull names),
Delphi is not case sensitive and this create a problem on a
computation.

If you are interessed I will send you my simple Pascal Code.
thank you again