Information och instruktioner nedan är hämtat från detta blogg-inlägg:
http://mathgis.blogspot.com/2008/01/det ... -data.htmlLadda hem program shp2text från följande adress :
http://www.obviously.com/gis/shp2text/shp2text.zipÖppna konsollen (cmd eller command från windows start-meny)
shp2text --gpx indiana.shp
titta på shape-filen för att se vilken information den innehåller.
shp2text --gpx indiana.shp 1 0 > output.xml
Detta sparar shape-filen som en XML-fil med data från Field 1 som namn och Field 0 som polygon-info.
T.ex.
shp2text --gpx indiana.shp
Field 0: Type=Integer, Title=`FIPS', Width=5, Decimals=0
Field 1: Type=String, Title=`COUNTY', Width=16, Decimals=0
Field 2: Type=String, Title=`DOQ_DISK', Width=16, Decimals=0
...
shp2text --gpx indiana.shp 1 0 > output.xml
I min kod nedan har jag tagit en massa postnummerområden runt stockholm.
Code:
dataxml = Import["E:\\GIS\\Sweden-SWEREF99\\sthlmpnr.xml"];
(*extract postal codes*)
pnr = Drop[Cases[dataxml, XMLElement["name", __, {w_}] :> w, Infinity], 1];
(*extract {lat,lon} for each county boundary*)
wholes = Cases[dataxml, XMLElement["rte", __], Infinity];
counties = (Cases[#, XMLElement["rtept", __], Infinity] &) /@ wholes;
dataLatLon =
counties /. (XMLElement["rtept",
latLonRules_, {}] :> {{"lat", "lon"} /. latLonRules});
(*tidy up a little,of course,you can combine them into one line code*)
test = Flatten[#, 1] & /@ dataLatLon;
test = Map[ToExpression, test];
test = Map[Reverse, test, 2];
(*then draw it*)
Graphics[{EdgeForm[Gray], FaceForm[LightGreen], Polygon[test]}]
Code:
extrude[pts_, h_] :=
Module[{vb, vt, len = Length[pts], nh}, If[! NumericQ[h], nh = 0., nh = N@h];
vb = Table[{pts[[i, 1]], pts[[i, 2]], 0}, {i, len}];
vt = Table[{pts[[i, 1]], pts[[i, 2]], nh}, {i, len}];
GraphicsComplex[
Join[vb, vt], {Polygon[Range[len]],
Polygon[Append[
Table[{i, i + 1, len + i + 1, len + i}, {i, len - 1}], {len, 1, len + 1,
2 len}]], Polygon[Range[len + 1, 2 len]]}]]
Length[test]
Out[2] := 103
Kolumn 1 i tabellen kommer att innehålla en polygon (visas med Graphics[] eller Graphics3D[])
Kolumn 2 i tabellen kommer att innehålla data över hur mycket varje område skall höjas upp i 3 D - bilden
Kolumn 3 i tabellen kommer att bestämma vilken färg varje polygon/område skall ha
Code:
mapdata = Table[{Polygon[test[[x]]], x, x/200}, {x, 1, Length[test]}];
(* Draw a 3D map *)
Graphics3D[{EdgeForm[
Darker[ColorData["DarkRainbow"][N@#[[3]]]]],
FaceForm[ColorData["DarkRainbow"][N@#[[3]]]],
extrude[#[[1, 1]], #[[2]]]} & /@ mapdata, Boxed -> False,
Lighting -> "Neutral", BoxRatios -> {1, 1, 0.2}, ImageSize -> 800]
