({{, (RefLink[GeoGridUnitDistance, paclet:ref/GeoGridUnitDistance])[proj,loc,α]
gives the actual geo distance corresponding to a unit distance on the geo grid obtained with projection proj, evaluated in the limit of small displacement from location loc in direction α.}})

OpenersNotesSection(Spacer1)Details and Options

The inverse of geo grid unit distance is also known as point scale or particular scale.

(RefLink[GeoGridUnitDistance, paclet:ref/GeoGridUnitDistance]) combines the global nominal scale (the shrinking factor of the geo model to the reference model of the map, traditionally denoted as 1:125000 and similar) and the local distortion of scale induced by the geo projection.

The result of (RefLink[GeoGridUnitDistance, paclet:ref/GeoGridUnitDistance])[…] corresponds to a ratio of (RefLink[Quantity, paclet:ref/Quantity]) geo distances on the geo model (of Earth or any other body) and dimensionless distances of the projected geo grid.

Points of maps including large regions ("large-scale maps") correspond to large values of geo distance scale, and points of "small-scale maps" correspond to small values of geo distance scale.

If geo grid unit distance is independent of azimuth at a point, then it is said to be isotropic at that point. A geo projection is conformal if and only if geo grid unit distance is isotropic at all points, though its value may still vary from point to point.

A geo projection can be given as a named projection "proj" with default parameters or as {"proj",params}, where "proj" is any of the entities of (RefLink[GeoProjectionData, paclet:ref/GeoProjectionData]) and params are parameter rules like "StandardParallels"->{33,60}. (RefLink[GeoProjectionData, paclet:ref/GeoProjectionData])["proj"] gives the default values of the parameters for the projection "proj".

The location loc can be given as a coordinate pair {lat,lon} in degrees, a geo position object like (RefLink[GeoPosition, paclet:ref/GeoPosition])[…] or (RefLink[GeoGridPosition, paclet:ref/GeoGridPosition])[…] or as a geo entity (RefLink[Entity, paclet:ref/Entity])[…].

The bearing or azimuthal direction α is an angle measured clockwise from true north. It can be given as a (RefLink[Quantity, paclet:ref/Quantity]) angle, as a number in degrees or as a named compass direction like "North", "NE" or "NEbE".

(RefLink[GeoGridUnitDistance, paclet:ref/GeoGridUnitDistance]) threads over its location and direction arguments.

Possible options of (RefLink[GeoGridUnitDistance, paclet:ref/GeoGridUnitDistance]) include:

OpenersPrimaryExamplesSection(Spacer1)Examples  (26)

OpenersExampleSection0(Spacer1)Basic Examples  (2)

Compute the unit geo distance induced by the Mercator projection at Copenhagen in the northeast direction:

scale=GeoGridUnitDistance["Mercator",(Entity["Copenhagen", Entity["City",{"Copenhagen","Copenhagen","Denmark"}], "Entity[\"City\", {\"Copenhagen\", \"Copenhagen\", \"Denmark\"}]", "city"]),"NE"]

Compute the final position of a geodesic of that length starting from Copenhagen in the northeast direction:

GeoDestination[(Entity["Copenhagen", Entity["City",{"Copenhagen","Copenhagen","Denmark"}], "Entity[\"City\", {\"Copenhagen\", \"Copenhagen\", \"Denmark\"}]", "city"]),{scale,45}]

Compare the geo path joining those locations with the map units and the scale bar:

GeoGraphics[{Red,Thick,GeoPath[{(Entity["Copenhagen", Entity["City",{"Copenhagen","Copenhagen","Denmark"}], "Entity[\"City\", {\"Copenhagen\", \"Copenhagen\", \"Denmark\"}]", "city"]),%}]},GeoProjection->"Mercator",GeoScaleBar->$UnitSystem,GeoRangePadding->(Quantity[30, "km", kilometers, "Kilometers"]),Frame->True]

If map units are made to correspond to inches, the traditional scale notation is 1:2468243, as given by the following:

DecimalForm[scale/Quantity[1,"Inches"]]

Use the Mollweide projection to construct a flat map based on a reference sphere of radius 6371:

proj={"Mollweide","ReferenceModel"->6371,"Centering"->{0,0}}

Then a unit of projected distance at London corresponds to a geo distance between these two values:

GeoGridUnitDistance[proj,(Entity["Paris", Entity["City",{"Paris","IleDeFrance","France"}], "Entity[\"City\", {\"Paris\", \"IleDeFrance\", \"France\"}]", "city"]),MinMax]

GeoGraphics[(Entity["Paris", Entity["City",{"Paris","IleDeFrance","France"}], "Entity[\"City\", {\"Paris\", \"IleDeFrance\", \"France\"}]", "city"]),GeoProjection->proj,GeoRange->Max[%],GeoScaleBar->$UnitSystem,Frame->True]

OpenersExampleSection1(Spacer1)Scope  (9)

Compute the geo grid unit distance for a geo projection at your current geo location in the northward direction:

GeoGridUnitDistance["Mollweide",Here,0]

These are the default values of the parameters of the "Mollweide" projection:

GeoProjectionData["Mollweide"]

Specify other values for the parameters of the projection:

GeoGridUnitDistance[{"Mollweide","ReferenceModel"->1000},Here,0]

GeoGridUnitDistance[{"Mollweide","Centering"->{-40,150}},Here,0]

Specify a location using a pair {lat,lon} in degrees:

GeoGridUnitDistance["Mercator",{-80.5,30.32},60]

Use locations with geo position heads:

GeoGridUnitDistance["Mercator",GeoPosition[{-80.5,30.32}],60]

GeoGridUnitDistance["Mercator",GeoPositionXYZ[{911681,533170,-6268998}],60]

GeoGridUnitDistance["Mercator",GeoGridPosition[{0.714285,-0.753047},"Albers"],60]

Specify a location using a geo (RefLink[Entity, paclet:ref/Entity]) object:

GeoGridUnitDistance["Albers",(Entity["Toronto", Entity["City",{"Toronto","Ontario","Canada"}], "Entity[\"City\", {\"Toronto\", \"Ontario\", \"Canada\"}]", "city"]),-30]

Compute the geo distance scale for a list of locations, all along the same direction:

GeoGridUnitDistance["Mollweide",{GeoPosition[{10,40}],GeoGridPosition[{-30.3,20},"Mercator"],{85,178}},35]

Convert the (RefLink[QuantityArray, paclet:ref/QuantityArray]) output into its normal form:

%//Normal

Specify the azimuth as a number of degrees:

GeoGridUnitDistance["Equirectangular",{80,-40},95]

Specify the same azimuth as a (RefLink[Quantity, paclet:ref/Quantity]) angle:

GeoGridUnitDistance["Equirectangular",{80,-40},Quantity[95,"AngularDegrees"]]

Use any other angular unit:

GeoGridUnitDistance["Equirectangular",{80,-40},Quantity[1,"Radians"]]

Compute the geo grid unit distance for a list of different azimuths at the same location:

GeoGridUnitDistance["Equirectangular",{80,-40},Range[0,90,15]]//Normal

The input can also be given as a (RefLink[QuantityArray, paclet:ref/QuantityArray]) object:

%==GeoGridUnitDistance["Equirectangular",{80,40},QuantityArray[Range[0,90,15],"AngularDegrees"]]

Compute the range of possible values of geo grid unit distance at a given location:

GeoGridUnitDistance["UrmayevI",(Entity["Moscow", Entity["City",{"Moscow","Moscow","Russia"}], "Entity[\"City\", {\"Moscow\", \"Moscow\", \"Russia\"}]", "city"]),MinMax]

Compare with the (RefLink[MinMax, paclet:ref/MinMax]) of values for each integer degree azimuth:

MinMax@GeoGridUnitDistance["UrmayevI",(Entity["Moscow", Entity["City",{"Moscow","Moscow","Russia"}], "Entity[\"City\", {\"Moscow\", \"Moscow\", \"Russia\"}]", "city"]),Range[-180,180]]

(RefLink[GeoGridUnitDistance, paclet:ref/GeoGridUnitDistance]) can efficiently process values for large numbers of locations:

locs=RandomGeoPosition["World",10^6];

azims=RandomReal[{-180,180},10^6];

GeoGridUnitDistance["Mollweide",locs,azims]//AbsoluteTiming

Select the same reference model and geo model to eliminate the effect of nominal scale:

scale[proj_,loc_,α_]:=QuantityMagnitude[GeoGridUnitDistance[{proj,"ReferenceModel"->1},loc,α,GeoModel->1],"Meters"]

The inverses of geo grid unit distance along meridians and parallels are traditionally denoted as h and k:

h[proj_,loc_]:=1/scale[proj,loc,0];
k[proj_,loc_]:=1/scale[proj,loc,90];

The behavior for these cylindrical projections is identical along parallels, but different along meridians:

Plot[{k["Mercator",{lat,0}],k["Equirectangular",{lat,0}],k["CylindricalEqualArea",{lat,0}]},{lat,0,90},PlotLegends->{"Mercator","Equirectangular","CylindricalEqualArea"}]

Plot[{h["Mercator",{lat,0}],h["Equirectangular",{lat,0}],h["CylindricalEqualArea",{lat,0}]},{lat,0,90},PlotLegends->{"Mercator","Equirectangular","CylindricalEqualArea"}]

OpenersExampleSection2(Spacer1)Options  (2)

OpenersExampleSubsection0(Spacer1)GeoModel  (1)

By default, (RefLink[GeoGridUnitDistance, paclet:ref/GeoGridUnitDistance]) returns values for Earth:

GeoGridUnitDistance["Bonne",{50,-100},90]

Performing the same computation on the corresponding point of the Moon returns smaller scales:

GeoGridUnitDistance["Bonne",{50,-100},90,GeoModel->"Moon"]

Choose a spherical model of specific radius:

GeoGridUnitDistance["Bonne",{50,-100},90,GeoModel->Quantity[2000,"Kilometers"]]

OpenersExampleSubsection1(Spacer1)UnitSystem  (1)

Select the unit system used to return the distance scale:

GeoGridUnitDistance["Wiechel",GeoPosition[{10,120}],20,UnitSystem->"Metric"]

GeoGridUnitDistance["Wiechel",GeoPosition[{10,120}],20,UnitSystem->"Imperial"]

They are the same value, but in different units:

%==%

OpenersExampleSection3(Spacer1)Properties & Relations  (11)

Take the Mercator projection on the default ellipsoidal model of Earth, a location and a direction:

proj={"Mercator","ReferenceModel"->"ITRF00"};
p=GeoPosition[{89,0}];
α=50;

Geo grid unit distance at p in direction α is the limit of the quotient of true and projected distances between p and a nearby point in direction α:

projectedDistance[p_,q_,proj_]:=EuclideanDistance[GeoGridPosition[p,proj]["GridXY"],GeoGridPosition[q,proj]["GridXY"]]

Table[dist/projectedDistance[p,GeoDestination[p,{dist,α}],proj],{dist,Quantity[PowerRange[10^5,0.01,0.1],"Meters"]}]

ListPlot[%,Joined->True,PlotMarkers->Automatic,Ticks->{None,Automatic}]

Compare with the computed value:

GeoGridUnitDistance[proj,p,α,GeoModel->"ITRF00",UnitSystem->"Metric"]

(RefLink[GeoGridUnitDistance, paclet:ref/GeoGridUnitDistance]) is periodic in azimuth with a period of 180 degrees:

proj={"Bonne","Centering"->{40,0}};
p=GeoPosition[{80,40}];

F[θ_Real]:=QuantityMagnitude[GeoGridUnitDistance[proj,p,θ],"Meters"]

Plot[F[θ],{θ,-180,360}]

Find the positions of one minimum and one maximum:

{FindMinimum[F[θ],θ],FindMaximum[F[θ],θ]}

{mindir,maxdir}=θ/.%[[All,2]]

Those correspond to the semiaxes of this Tissot ellipse:

GeoGraphics[{GeoCircle[p,1000],Red,Arrow@GeoPath[{p,1000,mindir}],Blue,Arrow@GeoPath[{p,1000,maxdir}]},GeoProjection->proj]

Geo distance scale can vary strongly with azimuth at a given point:

proj="Equirectangular";
p=(Entity["Svalbard", Entity["Country","Svalbard"], "Entity[\"Country\", \"Svalbard\"]", "country"])["CapitalCity"]

Plot[GeoGridUnitDistance[proj,p,θ],{θ,0,360}]

These are the minimum and maximum values attained:

GeoGridUnitDistance[proj,p,MinMax]

Geo grid unit distance can vary strongly from point to point for the same projection and azimuth:

locs=RandomGeoPosition["World",10000]

GeoGridUnitDistance["LambertAzimuthal",locs,20]

The result varies by more than two orders of magnitude:

MinMax[%]

Histogram[%%]

Geo distance scale is proportional to the geo model parameter:

GeoGridUnitDistance["Albers",GeoPosition[{80,30}],30,GeoModel->#]&/@Quantity[{1,10,100,1000,10000},"Kilometers"]

Geo grid unit distance is inversely proportional to the reference model and central scale parameters:

GeoGridUnitDistance[{"Albers","ReferenceModel"->#},GeoPosition[{80,30}],30]&/@{1,10,100,1000,10000}

GeoGridUnitDistance[{"Albers","CentralScaleFactor"->#},GeoPosition[{80,30}],30]&/@{1,10,100,1000,10000}

For an ellipsoidal projection, geo grid unit distance depends only slightly on the choice of datum or ellipsoid:

GeoGridUnitDistance[{"Mercator","ReferenceModel"->#},GeoPosition[{80,30}],30,GeoModel->#]&/@GeodesyData[]//MinMax

Equidistant projections have constant geo grid unit distance along special paths on the map:

projs=GeoProjectionData["Equidistant"]

For conic and cylindrical projections, this usually happens along meridians, at any location:

points=RandomGeoPosition["World",100]

Table[GeoGridUnitDistance[{proj,"ReferenceModel"->1},points,0]//MinMax,{proj,{"ConicEquidistant","CylindricalEquidistant","Equirectangular","GallIsographic"}}]

For the azimuthal equidistant projection, this happens for all directions from its centering:

GeoGridUnitDistance[{"AzimuthalEquidistant","Centering"->Here},Here,RandomReal[{-180,180},100]]//MinMax

For short distances, (RefLink[GeoDistance, paclet:ref/GeoDistance]) can be approximated as a product of projected distance by geo grid unit distance:

p=GeoPosition[{15,20}];
q=GeoPosition[{15.05,20.01}];

d=GeoDistance[p,q]

Compute projected distance in a given projection:

proj="Albers";

EuclideanDistance[GeoGridPosition[p,proj]["GridXY"],GeoGridPosition[q,proj]["GridXY"]]

Multiply by geo grid unit distance in the direction from p to q:

%GeoGridUnitDistance[proj,p,GeoDirection[p,q]]

The difference with the true result is smaller than 1%:

%/d

Compute geo distance along a meridian using any projection:

p=GeoPosition[{15,120}];
q=GeoPosition[{88,120}];

GeoGraphics[{Red,Thick,GeoPath[{p,q}]},GeoProjection->"VanDerGrintenII",Frame->True,GeoGridLines->(QuantityPostfix[10, "°", degrees, "AngularDegrees"]),GeoRangePadding->(QuantityPostfix[15, "°", degrees, "AngularDegrees"])]

Extract the projection selected by (RefLink[GeoGraphics, paclet:ref/GeoGraphics]) and compute the projected points:

proj=GeoProjection/.Options[%,GeoProjection]

gridp=GeoGridPosition[p,proj];gridq=GeoGridPosition[q,proj];

Here is the geo grid unit distance along the meridian, as a function of the projected y coordinate:

f[y_Real]:=GeoGridUnitDistance[proj,GeoGridPosition[{0,y},proj],0,GeoModel->"ITRF00"]

Compute the distance through a numerical integration:

NIntegrate[f[y],{y,gridp["GridY"],gridq["GridY"]}]

Compare with the unprojected geo distance:

%-GeoDistance[{p,q}]

Compare intervals of geo grid unit distance for different projections at the same point:

projs={"EckertGreifendorff","Bonne","Albers","Mollweide","Equirectangular"};

p=GeoPosition[{15,20}];

NumberLinePlot[QuantityMagnitude[GeoGridUnitDistance[{#,"ReferenceModel"->1},p,Interval],"Kilometers"]&/@projs,PlotLegends->projs]

The geo grid unit distance in conformal projections, like Mercator, is isotropic (does not depend on azimuth):

sphconformal={ "AugustEpicycloidal","EquatorialStereographic","HerschelConicConformal","Lagrange","LambertConicConformal","Littrow","Mercator","ObliqueMercator","PeirceQuincuncial","Stereographic","TransverseMercator"};

GeoProjectionData[#,"Conformal"]&/@sphconformal

GeoGridUnitDistance[#,GeoPosition[(Entity["New York City", Entity["City",{"NewYork","NewYork","UnitedStates"}], "Entity[\"City\", {\"NewYork\", \"NewYork\", \"UnitedStates\"}]", "city"])],MinMax]&/@sphconformal

The actual value of the projected unit distance varies from point to point for any given projection:

GeoGridUnitDistance["Mercator",{80,0},MinMax]

GeoGridUnitDistance["Mercator",{40,50},MinMax]

GeoGridUnitDistance["Mercator",{70,-150},MinMax]

Both isotropy and the dependence on latitude are clear in a map showing Tissot indicatrices:

GeoGraphics[Table[GeoDisk[{lat,lon},Quantity[300,"Miles"]],{lat, -80,80,20},{lon,-150,150,30}],GeoProjection->"Mercator"]

OpenersExampleSection4(Spacer1)Possible Issues  (1)

If a geo location cannot be projected, then the geo grid unit distance cannot be computed either:

GeoGridPosition[GeoPosition[{20,120}],"Orthographic"]

GeoGridUnitDistance["Orthographic",GeoPosition[{20,120}],90]

This location is not on the half-Earth covered by the "Orthographic" projection with default center:

GeoGraphics[GeoRange->"World",GeoProjection->"Orthographic",GeoGridLines->Automatic]

OpenersExampleSection5(Spacer1)Neat Examples  (1)

Compare the geo grid unit distance in the Mercator projection at different latitudes:

ratio[lat_]:=QuantityMagnitude[GeoGridUnitDistance["Mercator",{lat,0},90,GeoModel->180/Pi],"Meters"];

p[lat_,lon_]:={ratio[lat]lon,lat};

points=Table[p[lat,lon],{lat,0,88,11},{lon,-180,180,30}];

The following diagram shows at the bottom a scale at lower latitudes and how it changes when it is projected at higher latitudes:

Graphics[{Point[Flatten[points,1]],Line/@points,Thickness[.002],Red,Line/@Transpose[points]},Frame->True]

Sometimes this diagram is presented using the inverse ratio:

ip[lat_,lon_]:={1/ratio[lat]lon,lat};

ipoints=Table[ip[lat,lon],{lat,0,77,11},{lon,-180,180,30}];

Graphics[{Point[Flatten[ipoints,1]],Line/@ipoints,Thickness[.002],Red,Line/@Transpose[ipoints]},Frame->True,AspectRatio->1/5]

({{DocumentationFooterIconsRelatedFunction.png, ({{See Also}, {(RefLink[GeoDistance, paclet:ref/GeoDistance])  ▪  (RefLink[GeoGridUnitArea, paclet:ref/GeoGridUnitArea])  ▪  (RefLink[GeoGridDirectionDifference, paclet:ref/GeoGridDirectionDifference])  ▪  (RefLink[GeoProjectionData, paclet:ref/GeoProjectionData])  ▪  (RefLink[GeoScaleBar, paclet:ref/GeoScaleBar])}})}})

({{DocumentationFooterIconsRelatedGuide.png, ({{Related Guides}, {▪(RefLinkPlain[Geodesy, paclet:guide/Geodesy])}, {▪(RefLinkPlain[Maps & Cartography, paclet:guide/MapsAndCartography])}})}})