. Introduction

One of the main sources of cadastral activities, engineering projects, geographical information system (GIS) applications to be accurate, reliable and sustainable throughout the country is the unique projection coordinates used in scaled maps. In this context, it has been aimed to reduce the distortions occurring in maps to minimum levels by using correction formulas. This is achieved by scale and map harmony. The fact that the coordinate system is not single causes more than one start of coordinates throughout the country. In the zone system, only two neighbouring zones can be converted. If there are more than two zones, it is necessary to switch to another projection system. This leads to an increase in processing load and loss of time due to user errors.

Factors such as geodetic location, size and shape have been taken into consideration in order to minimise distortions in the projection selection of countries. However, due to both military and political aims, projections’ choice for these criteria has been taken into account and the Universal Transverse Mercator (UTM) system is preferred in Turkey. Turkey is a country that extends in the east–west direction. Therefore, the preference of UTM with increased projection distortion in this direction caused the use of more than one UTM zone in the country. Therefore, coordinate unity cannot be established in the ongoing Turkey's National Geographic Information System, the National Spatial and Spatial Data Infrastructure (NSDI) studies and engineering projects exceeding one slice. The aim of this study is to determine medium- and large-scale maps in the country and a projection with a single non-UTM coordinate origin for these projects. Therefore, the most suitable projection system for Turkey is the Lambert conformal conical (LCC2) projection type, in which the distortion values increase in the north–south direction. The studies of the countries on the selection of different projections other than UTM are available in the literature (Bugayevskiy and Snyder, 1995; Cory et al., 2001; Dennis, 2018; Habib, 2008; Hartzell et al., 2002; Huryeu and Padshyvalau, 2007; Vaníček and Najafi-Alamdari, 2004; Veverka, 2004). In this article, it has been aimed to determine a projection system by examining the applicability of LCC2 and external projection systems for Turkey, which lies in the east–west direction, parallel to latitude and located in the mid–latitude belt between the equator and the pole. Thus, a single-beginning and single-zone projection system throughout the country, along both the north–south and east–west directions, and a single coordinate system in all medium- and large-scale sheets are defined.

. Materials and Methods

Usage and disadvantages of UTM in Turkey

In Turkey, base maps in the scale of 1:25,000 are produced in Gauss–Kruger (GK) mapping system with 6° zones. On the other hand, large-scale cadastral maps in the scale of 1:1000 and standard topo-graphic maps in the scale of 1:5000 and 1:10,000 are produced in GK mapping system with 3° zones (LSMMIPR, 2020).

For Turkey, while the central meridians of GK projection with 3°-wide zones (MUTM) are 27°, 30°, 33°, 36°, 39°, 42° and 45°, these central meridians are 27°, 33°, 39° and 45° for GK projection with 6°-wide zones (UTM). There are four zones for the UTM system and the numbers of these zones are 35, 36, 37 and 38. In order to avoid misunderstandings, zone numbers are added before the easting value in UTM coordinates (Figure 1).

Figure 1

UTM zones for Turkey


In engineering projects that do not fit into a single zone, if there are overlapping parts with the neighbour zones in their 0.5° or 1° overlap area, the coordinate system is obtained by doing a zone transformation. If the application area overflows the overlap area, zone transformation would not be a suitable process. In this situation, it is not proper to work with a single coordinate system. A transformation from the UTM to different coordinate systems has to be done. The transformation and switching to other coordinate systems cause errors and loss of time because the users do not possess sufficient information. Therefore, the UTM is not sufficient for engineering projects that do not fit into a single zone. Hence, different representation methods should be used for these projects. However, the methods that implement different approaches other than the classical zone width approach have taken their place in the literature (Yildirim, 2004).

Projection methods and data set for determination of distortions

Given Turkey's geodetic location and considering the origin of a single coordinate system, projection zone width is increased (l ± 3.5°) UTM, Composite (CP), LCC2 and Double (DP). 585 test points have been selected in GRS80 ellipsoid and ITRF96 datum reference surface by creating a geographical grid network (0.5° × 0.5°) within the geographical borders of Turkey (B = [35.5°N–42.5°N]; L = [26°E–45°E]). In order to distribute the distortions that occur in all projections evenly, the starting coordinate has been taken as the centre of the geodetic boundaries (B0 = 39°N; L0 = 35.5°E). Thus, it has been aimed to examine the distortions that may occur in the border regions as they move away from the beginning. Common formulas have been used for the arc to chord (Tt) and distance (Ss) distortions in each projection. For these distortions, each test point has been chosen as a starting point and 45° azimuth angle and 1-km edges have been produced from this point. The coordinates of the second point of these edges have been calculated with the help of the direct solution of geodetic (Vincenty, 1975). After the calculation of projection distortions, arc to chord (Tt) map and distance (Ss) map correction values have been calculated with the help of projection coordinates. Nomenclature of variables used in projection algorithms has been given in the list below:

(B, L)Ellipsoid geographical latitude and longitude;
(φ, λ)Sphere geographical latitude and longitude;
(B0, L0)Central meridian (latitude) and longitude of central meridian;
(x, y)Projected coordinates or plane coordinates;
(N, E)Northing, Easting map of scale or sheet coordinates;
qIsometric latitude for ellipsoid: atanh(sin B) – e atanh(e sin B);
ωIsometric latitude for sphere: atanh(sin φ);
q, lIsometric coordinates;
AEllipsoid azimuth;
SEllipsoid geodetic or distance on ellipsoid;
sPlane or projection distance;
tPlane or projection azimuth;
γGrid convergence of ellipsoid;
(Tt)The arc to chord distortion of projection (At–γ);
(Ss)Distance distortion of projection (Ss);
(Tt)mapThe arc to chord correction calculated from (x, y);
(Ss)mapDistance correction calculated from (x, y);
aSemi-major axis of ellipsoid;
bSemi-minor axis of ellipsoid;
cPolar radius of curvature;
e2Eccentricity of ellipsoid squared;
e′2Second eccentricity of ellipsoid squared;
e2.71828182845904523. . . Euler's number;
GMeridian arc length from equator to latitude, meridian distance;
MRadius of curvature in the meridian:
c/(1 + e2 cos2 B)3/2;
NRadius of curvature in the prime vertical: c/(1 + e2 cos2 B)1/2;
RRadius of Gauss sphere: (MN)1/2 = c/(1 + e2 cos2 B);
rRadius of curvature of the parallel (N cos B);
η2e2 cos2 B;
ttan B;
mPoint grid scale factor;
m0Grid scale factor assigned to central meridian (longitude) or parallels (latitude);
(B1, B2)Standard parallels for Lambert conformal conic.

UTM methods with increased zone width

Inspired by the Mercator projection that had been developed in the 16th century, J. H. Lambert improved the Transversal Mercator (TM) projection for sphere in the 17th century. TM has been developed by C. F. Gauss in the early 19th century, but published by Schreiber and Kruger for ellipsoid. Therefore, TM is also used as GK. Later, ellipsoid formulas were re-examined by Hristow and the univariate and bivariate power series used today were obtained. The current usage type (Easting and Northing) of UTM was developed by the US Army in 1947. Many countries worldwide use the UTM system in mapping applications. Due to the increase in distortion in the east–west direction of the UTM system, 6° zone meridians have been used. In addition, UTM can be modified for large-scale maps and 3° zone meridians can be used regionally (Deakin et al., 2010; Snyder, 1987).

The formulas used in UTM could only be used for the maximum l = ±3.5° from the initial longitude chosen in the 19th century due to the lack of computational techniques. However, with the development of mathematical computer software for numerical computing (MATLAB, Maxima, Mathematica, Maple, etc.), new solutions (Bermejo-Solera and Otero, 2009; Bowring, 1990, 1993a,b; Deakin et al., 2010; Engsager and Poder, 2007; Grafarend, 1995; Guo et al., 2020; Ingwersen, 1996; Karney, 2011; Klotz, 1993; Lee, 1962, 1963; Mittermayer, 1993; Thompson, 1975; Turiño, 2004, 2008) have been developed for l > ±3.5° by increasing the number of coefficients in UTM. The coefficients used in the calculations are presented ready to the user. Mittermayer's method has been preferred in this study as it provides ease of calculation. The correlations in this method are given below:

where G90 symbolises the meridian arc length between 0° and 90° latitudes. Meridian convergence and differential scale factor used in determining distortions are calculated with the help of the following equations:

The bk coefficients in the Equation (3) and (10), are given in Table 1, in the range of 0–9. The geodetic boundaries of convergence of distortion and meridional in Turkey obtained using this method are shown in Figure 2.

Table 1

bk coefficients

kHayford (ED50 Datum)GRS80 (ITRF Datum)
Figure 2

Distortion sizes for Mittermayer method: (a) scale factor, counter interval 0.0005; (b) meridian convergence, contour intervals 1°; (c) the arc to chord (Tt) distortion, contour intervals 0.2°; (d) the distance (Ss) distortion, contour intervals 75 cm.


Lambert conformal conic projection

J. H. Lambert developed a single and double standard parallel conical projection for both sphere and ellipsoid in the 17th century. Double standard parallel (B1; B2) LCC2 projection is preferred to work with a single coordinate system in the entire country in order to reduce the distortions that occur from the selected B0 towards north and south in the x-axis direction. Projection coordinate data are calculated with the help of the following equations:


The meridian convergence and differential scale factor used in determining the distortions are calculated with the help of the following equations (Bugayevskiy and Snyder, 1995; Grossmann, 1976; Hooijberg, 2012; Thomas, 1952; Yang et al., 1999):


The selection of standard parallels varies according to the shape and size of the study area to be projected. Therefore, the developed equations are as follows:

where Bs and Bn symbolise the southern and northern latitude boundaries of the study area respectively. Also, k is a coefficient that changes according to the shape of the region. If the region is rectangular, k = 5; if it is a circle, k = 4; and if it is close to a rhombus, k = 3 (Bugayevskiy and Snyder, 1995). Standard parallels B1 = 37°30′N; B2 = 40°30′N and k = 5 are taken for Turkey. Calculated distortions and meridian convergence are shown in Figure 3.

Figure 3

Distortion sizes for LCC2: (a) scale factor, counter interval 0.0005; (b) meridian convergence, contour intervals 1°; (c) the arc to chord (Tt) distortion, contour intervals 0.2″; (d) the distance (Ss) distortion, contour intervals 40 cm.


Composite projection (CP)

The CP transcribes the map's geometry to scale, to the map's height-to-width ratio and to the central latitude of the shown area by replacing projections and adjusting their parameters. The CP indicates the entire globe including poles; it portrays continents or larger countries with less distortion and it can morph to the web Mercator projection for maps displaying small areas. By mixing different projections in order to reduce the deformation in scaled maps, applications have been developed for small-scale map applications (Jenny, 2012; Jenny and Šavrič, 2018; Šavrič and Jenny, 2014) as well as large- and medium-scale maps.

CPs (Gojamanov and Ismayilov, 2019; Huryeu and Padshyvalau, 2007, 2008; Padshyvalau et al., 2005; Pędzich, 2005) including Lagrange polyconic projection (Bugayevskiy and Snyder, 1995; Yang et al., 1999) and Chebyshev–Grave criteria (Milnor, 1969; Nestorov, 1997) are defined. CP is a type of projection that minimises distortions that increase in all directions as it moves away from the selected coordinate origin and aims to provide users with formula and calculation convenience. Basic formulas for CP are given below:

where (x1,y1) and (x2,y2) values are the GK and LCC2 projection co-ordinates used in this study, respectively. Considering the physical size of Turkey's geodetic boundaries, it is seen that there are distortions in both north–south and east–west directions. These two projection methods have been chosen to reduce these distortions in both directions. k1 and k2 are mixed projection coefficients and calculated using the following equation:

The value of n in the Equation (23) symbolises the number of points (585 points) calculated in practice. Coordinates of the second projection can also be used in the calculation of dk. In this study, the dk coefficients have been tested in both projections and LCC2 has been preferred because it reduces distortions. In order to prefer CP, the projection equations (x1,y1) and (x2,y2) in Equation (6) should be calculated. For this, isometric coordinate parameter pairs can be written in terms of complex forms in conform projections (Gojamanov and Ismayilov, 2019; Huryeu and Padshyvalau, 2007, 2008; Padshyvalau et al., 2005; Pędzich, 2005):


The following equations are used to avoid calculation overhead in complex functions:


Since the GK mapping accepts the equator as the beginning of the x-axis, the G0 value must be calculated. On the other hand, in the LCC2 description, since the x-axis is accepted as the beginning of the B0 latitude, it is not necessary to calculate the G0 value (G0 = 0). The ak coefficients are obtained by using the Taylor series at latitude B0 selected by taking l = 0 in accordance with the conformal description conditions of the analytical function derivative (Grossmann, 1976):


The values of δ and R0¯ in the equations given above for LCC2 are calculated using the Equations (1117). The effect of the ak coefficients on the (x, y) coordinates falls below 1 mm after the 10th iteration, considering the geographical borders of Turkey. These coefficients are difficult to derive and calculate manually. MATLAB software has been used for this process. These coefficients have been calculated using a number of functions (syms, diff and complex) available in the MATLAB library. The ak coefficients for CP are calculated by using Equations (2629) or Equations (3032) instead of Equations (21) and (22) (Gojamanov and Ismayilov, 2019; Huryeu and Padshyvalau, 2007, 2008; Padshyvalau et al., 2005; Pędzich, 2005):


Thus, the ease of calculation for CP is presented to the users with the coefficients in Table 2.

Table 2

ak coefficients in metres


kHayford (ED50 Datum)GRS80 (ITRF Datum)Hayford (ED50 Datum)GRS80 (ITRF Datum)Hayford (ED50 Datum)GRS80 (ITRF Datum)

The meridian convergence and differential scale factor used in determination of the distortions are calculated as follows.


The distortions and meridian convergence calculated are shown in Figure 4.

Figure 4

Distortion sizes for CP: (a) scale factor, counter interval 0.0005; (b) meridian convergence, contour intervals 1°; (c) the arc to chord (Tt) distortion, contour intervals 0.2″; (d) the distance (Ss) distortion, contour intervals 40 cm.


Conical double projection from ellipsoid to sphere and sphere to plane (DP)

The conformal projection of ellipsoid to sphere is fixed and there are many algorithms for conformal projection from sphere to plane. DP consists of two stages: firstly, from the isometric coordinates in the ellipsoid (q, l) to the isometric coordinates in the sphere (ω,Δλ), and then from the geodetic coordinates in the sphere (φ,λ) to the Gaussian sphere (Grossmann, 1976).


The double standard parallel projection equations from sphere to plane are as follows:


The meridian convergence and differential scale factor used in determination of distortions are calculated as follows (Bugayevskiy and Snyder, 1995; Grossmann, 1976):


The calculated distortions and meridian convergence are shown in Figure 5.

Figure 5

Distortion sizes for double projection: (a) scale factor, counter interval 0.0005; (b) meridian convergence, contour intervals 1°; (c) the arc to chord (Tt) distortion, contour intervals 0.2″; (d) the distance (Ss) distortion, contour intervals 40 cm.


Examination of the projections and a new proposal for Turkey

Considering the distortions in Figures 25, it is seen that the distortion in the Mittermayer method increases as it moves away from the starting meridian in the east–west direction compared to other methods. This means distortion of the size of the map in the whole country in large- and medium-sized maps. When the deformations in the other three methods are examined, it is seen that the edge distortion of CP gives better results between standard parallels compared to LCC2 and DP. In this study, for the first time in Turkey, CPU, GK and LCC2 depictions are used together. The combination of GK and LCC2 equations for the CP depiction equations can be seen as an additional computational burden on users. However, considering the geographical borders of Turkey, the coefficients to be used in the solution of Equations (21, 22) and (2629) are calculated and presented to the users in Tables 2 and 3. Thus, both ease of calculation is provided to the users and the edge distortion between the latitudes (37°30′ −40°30′) falls below 20 cm. The coefficients given in Tables 3 and 4 for CP are valid only for B0 = 39°. For the use of other countries, it is necessary to calculate a new average latitude value of the study area. Calculations for some points for CP are given in Table 5.

Table 3

bk coefficients


kHayford (ED50 Datum)GRS80 (ITRF Datum)Hayford (ED50 Datum)GRS80 (ITRF Datum)Hayford (ED50 Datum)GRS80 (ITRF Datum)
Table 4

Coefficients for correction formulas

Table 5

Numerical values for CP (for 1 km) (k1 = 0.001, k2 = 1 − k1)

B [°]L [°]XCP[m]YCP[m]Mm0=1mγ [°](Tt) [″](T-t)mapTt[m](Ss) [m](S-s)mapSs [m]

Calculation of geodetic coordinates from CP inverse projection coordinates

It is necessary to know the geodetic coordinates for transformation from the map coordinates of the newly defined system to the old system and for transition to geocentric coordinates in the Global Navigation Sattellite System (GNSS) application. This is done by inverse transformation from either GK or LCC2 projection coordinates (Equations (4144)):


For the value of B0 = 39°, the bk coefficients are calculated using the formulas given below:


All calculated bk coefficients are given in Table 3. Calculation of the values of δ and R0¯ used in the above equations are given in Equations (1117). The following equations are used for transition from isometric latitude to geodetic latitude given in Equations (5558) (Kaya, 1994):


Also, bk coefficients can be obtained by using Equation (66):


The reverse conversion process is performed using Equations (5558) and (66). This process is more practical for users. The coefficients calculated from Equations (5963) and (66) are shown in Table 3.

Correction formulas of CP

In engineering projects where distance and direction are measured, it is necessary to reduce lengths and direction. Since map coordinates are used in the projects, giving the reduction formulas over these coordinates provides ease of calculation for the users. The reduction formulas commonly used for conformal descriptions (Draheim, 1953) can also be used for CP:


The Δx(x2x1) and Δy(y2y1) reductions are the difference between the map coordinate values of the edge to be calculated. In order to provide ease of calculation to the users, tk and sk coefficients in Equations (6668) are calculated in B0 = 39° latitude and ITRF96 datum are given in Table 4. These coefficient values that can be used for up to 2 km in length can be measured in the land of Turkey's geodetic boundaries.

These differences Tt = (Tt) – (Tt)map and Ss = (Ss) – (Ss)map are calculated on the sides of 1 and 2 km and are given in Figure 6. It is observed that the edge distortion (Ss) in Equations (6668) is achieved with an accuracy of 1 cm for 1 km and less than 2 cm for 2 km. In addition, the map distortion (Tt) in Equations (6668) provides an accuracy of less than 0.003″ for 1 km and 0.006″ for 2 km. In this context, Equations (6668) are sufficient for the location accuracy of the angle and distance measurements made in the detail acquisition of projects using large-scale maps and in the application to the terrain.

Figure 6

Correction values based on correction formulas for accuracy analysis in CP: (a) Ss difference for 1 km, counter interval 2 mm; (b) Ss difference for 2 km, counter interval 2 mm; (c) Tt difference for 1 km, contour intervals 0.001″; (d) Tt difference for 2 km, contour intervals 0.001″.


Choice of scale correction factor (m0) and false origin for Turkey

Since large-scale maps (1:1000–1:10,000) cover small areas, m0 is taken as 1 for the whole country (Figure 4a). Thus, in engineering projects that require detailed measurement from the field and application to the field, the land and the map remain the same in 1:1 scale. This provides the users with the ease of completing the application without any additional transformation process in the application of the data produced from large-scale maps to the terrain. However, this is not the case with medium-scale maps (1:25,000–1:250,000). In these maps, due to the area they occupy, distortion increases in the near regions of the depiction centre and as they move away from the depiction centre. For the CP representation, these maps also increase as distortions increase from the chosen initial latitude to the north and south, and particularly at the system edges. It may be possible to keep this growth within a certain boundary by preserving the conditions of representation, and especially to reduce it within the country borders. So, m could be slightly reduced for the whole country. For this, m0 must be equal to 1 at xmax/2. Thus, in the initial latitude, m0 will be slightly smaller than 1, but in the border latitudes (m − 1), the difference will be halved. xmax/2 corresponds to the latitude B = 40.50° and m0 = 1 for CP. The m0value can be chosen as 0.9995, the average value between B = 40.5° and Bmax = 42° in CP (Figure 7).

Figure 7

Scale correction factor for CP


In order to not work with negative coordinates in scaled maps, at the beginning of the coordinate of the description, negative x values south of latitude B0 = 39° and negative y values west of longitude L0 = 35.5° should be positive. For this, a fixed value (false) must be specified to add to all projection coordinates (xcp, ycp). Negative values of CP coordinates are xmax = −388 km and ymax = −861 km in Turkey (Table 5). Therefore, false northing = 500,000 m and false easting = 1,000,000 m values are determined for positive coordinates. Thus, scale maps for CP (scale map) or sheet (sheet) coordinates can be calculated in Turkey:


In these equations, m0 value is taken as 1 for large-scale maps and 0.9995 for medium- and small-scale maps.

The meanings of the parameters symbolised in Table 5 are explained in Section 2.2. This table has been produced by taking into account the latitude and longitude limit values of Turkey. Distortion values have been calculated in CP as a result of calculations made using MATLAB software.

. Results and Discussion

In this study, projections that can be used outside the UTM system with the beginning of seven meridian zones are examined in Turkey. As it is known in the UTM system, zone transformation is performed in projects in two zones, and transition to a different projection is performed for projects falling into more than two zones. This leads to an increase in working load and working errors caused by users. In this context, Mittermayer with increased meridian zone width, CP, double standard parallel LCC2 and ellipsoid to sphere, spherical to plane conformable DP, which are suitable for the geodetic location and shape of our country, except UTM, have been examined in detail.

In this study, the applicability of CP, which is used for small-scale maps among the discussed projections, for large- and medium-scale maps has been investigated for the first time. Also, for the first time in CP, GK and LCC2 projections have been mixed in order to reduce distortions in all directions. The distortions of the four projection methods have been calculated with the test points determined in the 0.5° × 0.5° geodetic grid covering the geodetic borders of Turkey. According to the results, it has been observed that CP especially has less length distortion compared to LCC2. Thus, it is proven that a projection may also use CP and LCC2, except for Turkey.

Some definitions need to be made in order to be able to use CP in country mapping system, medium- and large-scale maps and engineering projects, which include detail and application processes. These definitions are the transformation of geographic and projection coordinates, direction and length correction formulas, scale correction factor (m0), false easting and false northing, respectively. Also, it is defined as the coefficient in terms of providing convenience to users in this study in Turkey. It is sufficient to calculate these coefficients once according to the initial latitude (B0) and they are valid everywhere within the geodetic boundaries of the country. These coefficients are calculated in ITRF96 and ED50 datum and given to users in this study.

Geodetic–projection coordinate transformation in CP is calculated with the help of Equations (2629, 38). The coefficients are given in Table 2. In addition, the Geodetic–projection coordinate mapping is performed using Equations (5558, 66). Coefficients are given in Table 3. Direction and length corrections are calculated with Equations (6769) and the coefficients are given in Table 4. Equations (6769) are sufficient to determine the position accuracy of the angle and distance measurements made with the total station in places where GNSS cannot be measured in the detail acquisition and application of the projects produced on large-scale maps.

. Conclusions

As is the case in Turkey, slice transformation should be performed for projects in two slices in the UTM system and a different projection system should be used for projects with more than two slices. This situation causes the working load to increase and user-based errors to occur. In this study, projections that can be used outside of UTM are examined in Turkey. Mittermayer, mixed projection (CP), double standard parallel LCC2 and conformal DP to ellipsoidal sphere and from sphere to plane, suitable for Turkey's geographical location and shape, have been studied in detail.

According to the results obtained, map projection CP with the least distortion values in both east–west and north–south directions has been selected. With the CP selection, a single coordinate system has been determined for medium- and large-scale maps. Projection correction formulas, scale factor and false origin have been calculated for map coordinates in CP. These distortions are obtained with a difference of less than 1 cm for 1 km long sides and less than 0.003″ for the azimuth value of this side, when the correction formulas are used. In addition, in this study, equations have been produced in the form of coefficients in order to provide convenience to the users. These coefficients are calculated separately for ITRF and ED50 datums according to the initial latitude (B0) and it is sufficient to calculate once. These coefficients are valid everywhere within the geographical boundaries.