(* Location on a map *)

(* Reading blank map *)
(* If CounryData fails, then use ClearSystemCache[] *)
maintainGeoBlank := (
filecdata = ToFileName[{dirwork}, "countrydata.m"];
If[FileType[filecdata] === File, Get[filecdata]];
If[!ValueQ[cdQatar] || !ValueQ[cdBahrain] || !ValueQ[cdSaudiArabia] || !ValueQ[cdUnitedArabEmirates],
   cdQatar = CountryData["Qatar", "Coordinates"][[2]]; (* Sometimes crashing Mathematica *)
   cdBahrain = CountryData["Bahrain", "Coordinates"][[1]];
   cdSaudiArabia = CountryData["Saudi Arabia", "Coordinates"][[1]];
   cdUnitedArabEmirates = CountryData["United Arab Emirates", "Coordinates"][[2]];
   Save[filecdata, {cdQatar, cdBahrain, cdSaudiArabia, cdUnitedArabEmirates}];
  ];
toMercator[{lat_, lng_}] := {lng, 
   Log[Abs[Sec[lat*Degree] + Tan[lat*Degree]]]/Degree};
cs = cdQatar;
cs1 = toMercator /@ cs;
csa = cdBahrain;
csa1 = toMercator /@ csa;
csb = cdSaudiArabia;
csb1 = toMercator /@ csb;
csc = cdUnitedArabEmirates;
csc1 = toMercator /@ csc;
(*file = ToFileName[{dirwork}, "map-of-qatar400.jpg"];*)
file = ToFileName[{dirwork}, "map-of-qatar.jpg"];
map = Import[file];
(* Adjusting map to the border *)
{xc, yc} = Transpose[cs1];
{xmin, xmax} = {Min[xc], Max[xc]};
{ymin, ymax} = {Min[yc], Max[yc]};
{x0, y0} = {xmin + xmax, ymin + ymax}/2;
{dx, dy} = {xmax - xmin, ymax - ymin};
{xsize, ysize} = {400, 566};
(* Increasing size *)
dy = 1.18 dy;
dx = xsize dy/ysize;
(* shift *)
corr = {0.007, -0.015};
transf = ((# - {x0, y0})/{dx, dy} + {0.5, 0.5} - corr) &;
cs2 = transf /@ cs1;
Print["transf = ", Map[Evaluate, transf] // InputForm];
transf = {0.493 + 0.6740514291908127*(-51.18955 + #1), 0.515 + 0.4763614340571115*(-26.226486888492424 + #1)} &;
);

maintainGeoMap := (
If[Head[pts] =!= List || pts === {},
   pts = {{24.8931, 50.907}, {24.8989, 51.0134}, {24.75, 51.0763}, {24.7489, 
  51.0725}, {25.0421, 50.8401}, {25.0151, 51.0534}, {25.4672, 
  51.2126}, {25.0886, 51.1362}, {24.9629, 51.0064}, {25.3199, 
  51.5119}, {25.3187, 51.5094}, {25.3187, 51.5094}};
   Print["Warning: maintainGeoMap: using default pts (only for testing): ", pts]];
If[Head[fileexp] =!= String || fileexp === "",
   fileexp = "F:\\temp\\test.jpg";
   Print["Warning: maintainGeoMap: using default fileexp (only for testing): ", fileexp]];
(* Input: pts, fileexp *)
pts = toMercator /@ pts;
pts = {EdgeForm[Directive[Thickness[0.01], Red]], 
     RGBColor[1, 0, 1, 1], Disk[#, {0.025, 0.025}]} & /@ pts;
opts = Graphics[{Texture[map], 
    Polygon[cs1, VertexTextureCoordinates -> cs2], Thick, Line[cs1], 
    pts}] // AbsoluteOptions;
{xrange, yrange} = PlotRange /. opts;
{xmean, ymean} = Mean /@ {xrange, yrange};
{xrange, yrange} = {xmean, ymean} + 1.05 ({xrange, yrange} - {xmean, ymean});
graypol[csx1_] := {EdgeForm[{(*Thickness[0.006], *)GrayLevel[0.4]}], GrayLevel[0.6], Polygon[csx1]};
plt = Graphics[{Texture[map], graypol[csa1], graypol[csb1], 
   graypol[csc1], Polygon[cs1, VertexTextureCoordinates -> cs2], 
   Thick, Line[cs1], pts}, PlotRange -> {xrange, yrange}, Background -> GrayLevel[192./256]];
Export[fileexp, plt, ImageSize -> 200, "CompressionLevel" -> 0.1];
fileexp1 = StringReplace[fileexp, ".jpg" -> "-big.jpg"];
Export[fileexp1, plt, ImageSize -> 600, "CompressionLevel" -> 0.1];
);