### How To Compute & Use UTM Scale Factor?

Recently i do some work on computing point scale factor in utm coordinate system. The input of system is Easting,Northing,height and zone number in UTM (WGS84) and output is point scale factor.

First of all i should mention that the thing that we use as scale factor most of the time is Grid factor  (Combined Scale Factor that is the product of UTM Scale factor and Elevation Factor). I explain more.

Consider three coordinate systems.

System1:The system that we do observation in the surface of the earth with traditional surveying instrument (like total station and theodolites but not GPS)

System2: WGS84 (for example) or any other ellipsoid that made a coordinate system (we make a coordinate system using an ellipsoid by projecting the point into the ellipsoid plane and measure geodetic fi and landa angles and also height of the point)

System3: is a projection surface(UTM) . this is what made (fi,landa) to (x,y) and make it possible to use paper maps to represent the earth curved surface (Earth is a curved surface in all points)

for converting a length from system1 to susyem3 we should use two kind of scale factors.

for example you measure a length on earth with GPS (E,N,h) and want to convert it to local length . what should we do?

first you should compute a scale factor that bring the point of observation from elevation h  into the surface of ellipsoid (h=0). it is because all of computations is on the surface. the name of this kind of scale factor is "Elevation Factor" (EF). it's formula is a very simple relation.

EF=(Re/Re+h) where Re= mean earth radius =6371000m and "h" is ellipsoid height of point.

second scale factor is what that transform between projection coordinate system(UTM) and WGS84 coordinate system.it's formula (with negligible truncation) is:

k=k0(1+(landa-landa0)^2*cos(fi)^2/2) where:

k=scale factor

k0=central meridian scale factor =0.9996 (in UTM Coordinate system)

(fi,landa) are the position of the point in WGS84

landa0=longitude of central meridian in that zone in utm Coordinate System

and at last for computing "grid scale factor" or combined scale factor,  you should multiply k*EF=Grid Factor

Grid factor is used to convert length from surface of earth into projection system (UTM) and if you want convert UTM length into earth surface length you should divide length into the Grid Factor.

you are done here but if you want go further, in measuring a length in fact we have 2 points, Start point and End point. it is better that you compute point scale factor in Start and End (or Start Middle and End) and then use average (or Weighted average) Scale Factor. in short distances it is not problem to use ONE SF but in long distances you need to average.

I wrote a C# program that takes E.N,h and zone and compute (fi,landa) and then Scal factor,Elevation Factor and finally Grid factor (combined Factor). here is a source code.

For rebuilding project you should just make a form and put some textbox controls and give desired input:

using System;
using System.Collections.Generic;
using System.ComponentModel;
using System.Data;
using System.Drawing;
using System.Text;
using System.Windows.Forms;
using System.IO;
using netDxf;
using netDxf.Blocks;
using netDxf.Collections;
using netDxf.Entities;
using netDxf.Objects;
using netDxf.Tables;
using Group = netDxf.Objects.Group;
using Point = netDxf.Entities.Point;
using Attribute = netDxf.Entities.Attribute;
using Image = netDxf.Entities.Image;

namespace UTM
{

public partial class Form1 : Form
{
public int Cs = 1;//0==WGS84,1==UTM
public double eccSquared = 0.00669438; // eccentricity (0.081819191 ^ 2) WGS84
public double dEquatorialRadius = 6378137.0; // WGS84 (note above: varies from 6,356.750 km to 6,378.135 km)
public double dScaleFactor = 0.9996; // scale factor, used as k0
public double dDenominatorOfFlatteningRatio = 298.257223563;
public double
dCvtLiterHa2GalAc = 0.10691,   // l/ha -> gal/ac (dCvtLiter2Gal / dCvtHa2Ac)
dCvtLiterHa2OzAc = 13.68416,   // l/ha -> fl oz/ac (dCvtLiter2Oz / dCvtHa2Ac)
dCvtDeg2Rad = Math.PI / 180.0, // 0.0174532925199432957 ... remember: pi radians = 180 deg
dCvtRad2Deg = 180.0 / Math.PI; // 57.2957795130823208767 ...

# region UTM
public Form1()
{
InitializeComponent();
}
/// <summary>
/// figure the utm zone based solely on longitude (in degrees)
/// </summary>
/// <param name="dLon"></param>
/// <returns></returns>
public int iGetUtmZone(double dLat,
double dLon)
{
// Make sure the longitude is between -180.00 .. 179.9
// ... could force it as follows:
//     double LongTemp = (Long + 180) - ((int)((Long + 180) / 360)) * 360 - 180; // -180.00 .. 179.9;
if ((dLon < -180.0) || (dLon > 180.0))
return (0); // implies failure: 0 is an invalid utm zone

int iUTM_Zone_Num = (int)((180.0 + dLon) / 6.0 + 1.0);

if (dLat >= 56.0 && dLat < 64.0) // sw Norge (ie, zone 32V)
{
if (dLon >= 3.0 && dLon < 12.0)
iUTM_Zone_Num = 32;
}
else if (dLat >= 72.0 && dLat < 84.0) // special zones for Svalbard
{
if (dLon >= 0.0 && dLon < 9.0)
iUTM_Zone_Num = 31;
else if (dLon >= 9.0 && dLon < 21.0)
iUTM_Zone_Num = 33;
else if (dLon >= 21.0 && dLon < 33.0)
iUTM_Zone_Num = 35;
else if (dLon >= 33.0 && dLon < 42.0)
iUTM_Zone_Num = 37;
}

return (iUTM_Zone_Num);
}
/// <summary>
///
/// </summary>
/// <param name="Lat"></param>
/// <returns></returns>
private string sUtmLetterNS(double dLat)
{
if (dLat >= 0)
return ("N");
return ("S");
}
/// <summary>
///
/// </summary>
/// <param name="iZoneNum"></param>
/// <returns></returns>
public double dSet_CentralMeridian_from_UtmZone(int iZoneNum)
{
double dCentralMeridian = (183.0 - (6.0 * (double)iZoneNum)) * -1.0;
return (dCentralMeridian);
}
public int iLatLon2UTM(double dLat,
double dLon,
ref double UTMNorthing,
ref double UTMEasting,
ref string sUtmZone)
{
// if deviation from WGS84 is desired, do this (after loading the array, duh):
//  eccSquared = ellipsoid[index-of-desired-reference-ellipsoid].eccentricitySquared;

// use dLonWork to make sure the longitude is between -180.00 .. 179.9
double dLonWork = (dLon + 180) - ((int)((dLon + 180) / 360)) * 360 - 180; // -180.00 .. 179.9;

int iUTM_Zone_Num = this.iGetUtmZone(dLat, dLonWork);
// set the resultant UTM Zone string
sUtmZone = iUTM_Zone_Num.ToString() + sUtmLetterNS(dLat);

// set central meridian
double dCentralMeridian = dSet_CentralMeridian_from_UtmZone(iUTM_Zone_Num);

double eccPrimeSquared = eccSquared / (1 - eccSquared);

double M = dEquatorialRadius * ((1 - eccSquared / 4 - 3 * eccSquared * eccSquared / 64 - 5 * eccSquared * eccSquared * eccSquared / 256) * dLatRad
- (3 * eccSquared / 8 + 3 * eccSquared * eccSquared / 32 + 45 * eccSquared * eccSquared * eccSquared / 1024) * Math.Sin(2 * dLatRad)
+ (15 * eccSquared * eccSquared / 256 + 45 * eccSquared * eccSquared * eccSquared / 1024) * Math.Sin(4 * dLatRad)
- (35 * eccSquared * eccSquared * eccSquared / 3072) * Math.Sin(6 * dLatRad));

UTMEasting = (double)(dScaleFactor * N * (A + (1 - T + C) * A * A * A / 6
+ (5 - 18 * T + T * T + 72 * C - 58 * eccPrimeSquared) * A * A * A * A * A / 120)
+ 500000.0);

UTMNorthing = (double)(dScaleFactor * (M + N * Math.Tan(dLatRad) * (A * A / 2 + (5 - T + 9 * C + 4 * C * C) * A * A * A * A / 24
+ (61 - 58 * T + T * T + 600 * C - 330 * eccPrimeSquared) * A * A * A * A * A * A / 720)));
if (dLat < 0)
UTMNorthing += 10000000.0; // 10000000 meter offset for southern hemisphere

return (0);
}
public int iUTM2LatLon(double UTMNorthing,
double UTMEasting,
string sUTMZone,  // expected format "12N"
ref double dLat,
ref double dLon)
{
// if deviation from WGS84 is desired, do this (after loading the array, duh):
//  eccSquared = ellipsoid[index-of-desired-reference-ellipsoid].eccentricitySquared;

// populate North/South
char cZoneLetter = sUTMZone[sUTMZone.Length - 1];
bool bNorthernHemisphere = (cZoneLetter >= 'N');

string sZoneNum = sUTMZone.Substring(0, sUTMZone.Length - 1);
int iZoneNumber = Convert.ToInt32(sZoneNum);

double x = UTMEasting - 500000.0; //remove 500,000 meter offset for longitude
double y = UTMNorthing;
if (!bNorthernHemisphere) // point is in southern hemisphere
y -= 10000000.0; // remove 10,000,000 meter offset used for southern hemisphere

double dLongOrigin = (iZoneNumber - 1) * 6 - 180 + 3; // +3 puts origin in middle of zone

double eccPrimeSquared = (eccSquared) / (1 - eccSquared);

double M = y / dScaleFactor;
double mu = M / (dEquatorialRadius * (1 - eccSquared / 4 - 3 * eccSquared * eccSquared / 64 - 5 * eccSquared * eccSquared * eccSquared / 256));

double e1 = (1 - Math.Sqrt(1 - eccSquared)) / (1 + Math.Sqrt(1 - eccSquared));
double phi1Rad = mu + (3 * e1 / 2 - 27 * e1 * e1 * e1 / 32) * Math.Sin(2 * mu)
+ (21 * e1 * e1 / 16 - 55 * e1 * e1 * e1 * e1 / 32) * Math.Sin(4 * mu)
+ (151 * e1 * e1 * e1 / 96) * Math.Sin(6 * mu);
// convert to degrees

double R1 = dEquatorialRadius * (1 - eccSquared) / Math.Pow(1 - eccSquared * Math.Sin(phi1Rad) * Math.Sin(phi1Rad), 1.5);
double D = x / (N1 * dScaleFactor);

dLat = phi1Rad - (N1 * Math.Tan(phi1Rad) / R1) * (D * D / 2 - (5 + 3 * T1 + 10 * C1 - 4 * C1 * C1 - 9 * eccPrimeSquared) * D * D * D * D / 24
+ (61 + 90 * T1 + 298 * C1 + 45 * T1 * T1 - 252 * eccPrimeSquared - 3 * C1 * C1) * D * D * D * D * D * D / 720);
// convert to degrees

dLon = (D - (1 + 2 * T1 + C1) * D * D * D / 6 + (5 - 2 * C1 + 28 * T1 - 3 * C1 * C1 + 8 * eccPrimeSquared + 24 * T1 * T1)
* D * D * D * D * D / 120) / Math.Cos(phi1Rad);
// convert to degrees
dLon = dLongOrigin + dLon * dCvtRad2Deg;
return (0);
}
private void button1_Click(object sender, EventArgs e)
{
try
{
double dLat = 0;
double dLon = 0;
double UTMNorthing = 0;
double UTMEasting = 0;
string sUTMZone = "39N";

if (Cs == 1)
{

//cBaseCmnGIS icBaceCmnGIS = new cBaseCmnGIS();
UTMNorthing = System.Convert.ToDouble(txtN.Text);
UTMEasting = System.Convert.ToDouble(txtE.Text);
sUTMZone = txtZone.Text;  // expected format "12N"
int iZoneNum = Convert.ToInt32(sUTMZone.Substring(0, sUTMZone.Length - 1));

int iRes = iUTM2LatLon(UTMNorthing, UTMEasting, sUTMZone, ref dLat, ref dLon);

double d, m, s;
double sgn = 1;
if (dLon<0)
{
sgn = -1;
}
dLon *= sgn;
d = Math.Floor(dLat);
m = Math.Floor((dLat - d) * 60);
s = ((dLat - d) * 60 - m) * 60;//Math.Round(((dLat - d) * 60 - m) * 60, 4);
//lblfi.Text = "fi: " + d.ToString() + "°  " + m.ToString() + "'  " + s.ToString() + "\"";
txtFiD.Text = d.ToString();
txtFiM.Text = m.ToString();
txtFiS.Text = s.ToString();//Math.Round(s, 2).ToString();
d = Math.Floor(dLon);
m = Math.Floor((dLon - d) * 60);
s = ((dLon - d) * 60 - m) * 60;//Math.Round(((dLon - d) * 60 - m) * 60, 4);
//lblLanda.Text = "Landa: " + d.ToString() + "°  " + d.ToString()+ "'  " + s.ToString() + "\"";
txtLandaM.Text = m.ToString();
txtLandaS.Text = s.ToString(); //Math.Round(s, 3).ToString();
double k1, k2, k;//ScaleFactor
double dCentralMeridian = (183.0 - (6.0 * (double)iZoneNum)) * -1.0;
double d1 = dCvtDeg2Rad * (dCentralMeridian - dLon);
double d2 = Math.Cos(dLat * dCvtDeg2Rad);
k1 = dScaleFactor * (1 + d1 * d1 * d2 * d2 / 2);
k = k1 * k2;
txtProjectionSF.Text = k1.ToString();
txtEF.Text = k2.ToString();
txtCombinedSF.Text = k.ToString();
}
else if (Cs == 0)
{

if (string.IsNullOrEmpty(txtFiD.Text)) txtFiD.Text = "0";
if (string.IsNullOrEmpty(txtFiM.Text)) txtFiM.Text = "0";
if (string.IsNullOrEmpty(txtFiS.Text)) txtFiS.Text = "0";
if (string.IsNullOrEmpty(txtLandaM.Text)) txtLandaM.Text = "0";
if (string.IsNullOrEmpty(txtLandaS.Text)) txtLandaS.Text = "0";
dLat = System.Convert.ToDouble(txtFiD.Text) + System.Convert.ToDouble(txtFiM.Text) / 60 + System.Convert.ToDouble(txtFiS.Text) / 3600;
dLon = System.Convert.ToDouble(txtLandaD.Text) + System.Convert.ToDouble(txtLandaM.Text) / 60 + System.Convert.ToDouble(txtLandaS.Text) / 3600;
int Res = iLatLon2UTM(dLat, dLon, ref UTMNorthing, ref UTMEasting, ref sUTMZone);

txtN.Text = UTMNorthing.ToString();
txtE.Text = UTMEasting.ToString();
txtZone.Text = sUTMZone;

}
}
catch (Exception Err)
{
Console.WriteLine("{0} Exception caught.", Err);
}

}
private void Scalefactor() { }

/*  private void txtE_Enter(object sender, EventArgs e)
{
Cs = 1;
//MessageBox.Show("AA");
}

private void txtN_Enter(object sender, EventArgs e)
{
Cs = 1;
}

private void txth_Enter(object sender, EventArgs e)
{
Cs = 1;
}

private void txtZone_Enter(object sender, EventArgs e)
{
Cs = 1;
}

private void txtFiD_Enter(object sender, EventArgs e)
{
Cs = 0;
}

private void txtFiM_Enter(object sender, EventArgs e)
{
Cs = 0;
}

private void txtFiS_Enter(object sender, EventArgs e)
{
Cs = 0;
}

private void txtLandaD_Enter(object sender, EventArgs e)
{
Cs = 0;
}

private void txtLandaM_Enter(object sender, EventArgs e)
{
Cs = 0;
}

private void txtLandaS_Enter(object sender, EventArgs e)
{
Cs = 0;
}*/

private void gbUTM_Enter(object sender, EventArgs e)
{
Cs = 1;
}

private void gbWGS84_Enter(object sender, EventArgs e)
{
Cs = 0;
}

#endregion

#region DXF
/// <summary>
/// DXF Section
/// </summary>
string[] comboContent = { "col1", "col2", "col3", "col4", "col5", "col6", "Null" };
DxfDocument dxf;
string FileName;

private void Form1_Load(object sender, EventArgs e)
{
createTooltips();
setUpComboBoxes();

}

private void setUpComboBoxes()
{
try
{
// cmb_DXFVersion.DataSource = Enum.GetValues(typeof(DxfVersion));
string[] DXFVersions = Enum.GetNames(typeof(DxfVersion));
for (int i = 5; i <= 7; i++)

cmb_PtName.SelectedIndex = 0;
cmb_E.SelectedIndex = 1;
cmb_N.SelectedIndex = 2;
cmb_H.SelectedIndex = 3;
cmb_Code.SelectedIndex = comboContent.Length - 1;
cmb_Desc.SelectedIndex = comboContent.Length - 1;
cmb_DXFVersion.SelectedIndex = cmb_DXFVersion.Items.Count - 1;

}
catch
{
}
}

private void createTooltips()
{
try
{

// Create the ToolTip and associate with the Form container.
ToolTip toolTip1 = new ToolTip();

// Set up the delays for the ToolTip.
toolTip1.AutoPopDelay = 5000;
toolTip1.InitialDelay = 500;
toolTip1.ReshowDelay = 500;
// Force the ToolTip text to be displayed whether or not the form is active.
toolTip1.ShowAlways = true;

// Set up the ToolTip text for the Button and Checkbox.
toolTip1.SetToolTip(this.button1, "Compute UTM/LatLon and Scale Factors");
toolTip1.SetToolTip(this.btnBrowse, "Open *.txt|*.csv of Coordinates");
toolTip1.SetToolTip(this.txt_FilePath, "Shows the Current Coornate Cile Path");
toolTip1.SetToolTip(this.cmb_Deliminator, "Input Deliminator or select From List");
toolTip1.SetToolTip(this.cmb_PtName, "Select Column of Point Name");
toolTip1.SetToolTip(this.cmb_E, "Select Column of Easting of Points");
toolTip1.SetToolTip(this.cmb_N, "Select Column of Northing of Points");
toolTip1.SetToolTip(this.cmb_H, "Select Column of Height of Points");
toolTip1.SetToolTip(this.cmb_Code, "Select Column of Code of Points");
toolTip1.SetToolTip(this.cmb_Desc, "Select Column of Description of Points");
toolTip1.SetToolTip(this.txt_FontSize, "Sets Height of Text in DXF");
toolTip1.SetToolTip(this.cmb_DXFVersion, "Select DXF Version");
toolTip1.SetToolTip(this.chk2D, "Ignore H of Coordinate and insert All Object With H=0");
toolTip1.SetToolTip(this.chk_Binary, "If Checked a Binary DXF Will be Created");

toolTip1.SetToolTip(this.btn_Save, "Start to Create DXF");
toolTip1.SetToolTip(this.txt_SavePath, "Shows The Path of Created DXF");
toolTip1.SetToolTip(this.tabPage_DXF, "Gets a Coordinatex text File and Create a DXF File");
toolTip1.SetToolTip(this.tabPage_UTM, "Gets UTM/LatLon coordinate and Compute Scale Factor of Point");
}
catch
{
}
}

private void btnBrowse_Click(object sender, EventArgs e)
{
openFileDialog_Coordinate.FilterIndex = 1;

DialogResult result = openFileDialog_Coordinate.ShowDialog(); // Show the dialog.
if (result == DialogResult.OK) // Test result.
{
string file = openFileDialog_Coordinate.FileName;
txt_FilePath.Text = file;
FileName = openFileDialog_Coordinate.SafeFileName;

}

}

private void btn_Save_Click(object sender, EventArgs e)
{
try
{
dxf = new DxfDocument();
createObjexts();

string dxfFileName = FileName.Remove(FileName.Length - 3, 3) + "dxf";
string dxfPathName = txt_FilePath.Text.Remove(txt_FilePath.Text.Length - FileName.Length, FileName.Length) + dxfFileName;
// saving to file stream
using (FileStream fileStream = new FileStream(dxfPathName, FileMode.Create))
{
if (!dxf.Save(fileStream, chk_Binary.Checked))
{
throw new Exception("Error saving to file stream.");
}
txt_SavePath.Text = dxfPathName;
}
}
catch
{
}

}

private bool createObjexts()
{
txt_SavePath.Text = "";
try
{
char delimiter;
if (cmb_Deliminator.Text == "")
return false;
if (cmb_Deliminator.Text.Trim().Length == 1)
{
delimiter = Convert.ToChar(cmb_Deliminator.Text);
}
else if (cmb_Deliminator.SelectedIndex == 0)
{
delimiter = Convert.ToChar("\t");
}
else if (cmb_Deliminator.SelectedIndex == 1)
{
delimiter = Convert.ToChar(" ");
}

else
return false;
//char delimiter = Convert.ToChar(cmb_Deliminator.Text.Substring(0));
bool readPt_Name = !cmb_PtName.SelectedItem.ToString().ToUpper().Contains(comboContent[comboContent.Length - 1].ToUpper());
bool readCode = !cmb_Code.SelectedItem.ToString().ToUpper().Contains(comboContent[comboContent.Length - 1].ToUpper());
bool readDesc = !cmb_Desc.SelectedItem.ToString().ToUpper().Contains(comboContent[comboContent.Length - 1].ToUpper());
bool hasH = !cmb_H.SelectedItem.ToString().ToUpper().Contains(comboContent[comboContent.Length - 1].ToUpper());
double H = 0;

double txtOffset = Convert.ToDouble(txt_FontSize.Text) * 2;
short colorIndex;

foreach (string record in All_lines)
{
string[] record_split = record.Split(delimiter);
double E = Convert.ToDouble(record_split[cmb_E.SelectedIndex]);
double N = Convert.ToDouble(record_split[cmb_N.SelectedIndex]);

colorIndex = 70;
string coord_EN_layerName = "dxf_coord_EN";
addTexts(E.ToString(), E, N + txtOffset, H, coord_EN_layerName, TextAlignment.TopRight, colorIndex);
addTexts(N.ToString(), E, N, H, coord_EN_layerName, TextAlignment.MiddleRight, colorIndex += 10);

if (hasH)
{
string Htext = record_split[cmb_H.SelectedIndex];
if (chk2D.Checked == false)
H = Convert.ToDouble(Htext);
addTexts(Htext, E, N - txtOffset, H, "dxf_coord_EN", TextAlignment.BottomRight, colorIndex += 10);
}

string Name, Code, Desc;
colorIndex = 140;

{
Code = record_split[cmb_Code.SelectedIndex];
addTexts(Code, E, N + txtOffset, H, "dxf_Code", TextAlignment.TopLeft, colorIndex);
}
{
Name = record_split[cmb_PtName.SelectedIndex];
addTexts(Name, E, N, H, "dxf_PtName", TextAlignment.MiddleLeft, colorIndex += 10);
}
{
Desc = record_split[cmb_Desc.SelectedIndex];
addTexts(Desc, E, N - txtOffset, H, "dxf_Desc", TextAlignment.BottomLeft, colorIndex += 10);
}
}
}
catch
{
return false;
}
return true;
}

private bool addpoint(double E, double N, double H)
{
try
{
Point point = new Point(E, N, H);
point.Layer = new Layer("dxf_points");
//point.Color = new AciColor(30);
}
catch
{
return false;
}
return true;
}

private void addTexts(string txt, double E, double N, double H, string layerName, TextAlignment alignment, short colorIndex)
{
//text E
Vector3 vec3d = new Vector3(E, N, H);
//TextStyle style = new TextStyle("True type font", "Arial.ttf");
Text text = new Text(txt, vec3d, Convert.ToDouble(txt_FontSize.Text));
text.Layer = new Layer(layerName);
text.Layer.Color.Index = colorIndex;
text.Alignment = alignment;

}
#endregion

}
}

The precision is tested with Leica Geo Office tools (LGO) and i didn't see any problem.

In this little application you enter E,N,h in UTM coordinate system and also your zone. at the end of zone number you  should add N or S for North hemisphere or South hemisphere. then press compute.

The default location that are fill in fields (this is the location of my home ;-)) may be a good guide.

I thank many internet source code and articles.

If you have any question about similar problems i'll be so happy if i can help.