nst-guide / data Goto Github PK
View Code? Open in Web Editor NEWCreate map waypoints and layers from open data sources
Home Page: https://nst-guide.github.io/data/
License: GNU General Public License v3.0
Create map waypoints and layers from open data sources
Home Page: https://nst-guide.github.io/data/
License: GNU General Public License v3.0
There are parts of my recorded GPS track with bad data. For geotagging my photos, it would be good to
Project own track onto OSM track. This means that every photo will have a position on the trail.
Compare iPhone GPS coordinates to interpolated coords from watch
TODO
comment in 444ca1d. It's been assigned to @kylebarron because they committed the code.TODO
comment in 660e995. It's been assigned to @kylebarron because they committed the code.Current relation: https://www.openstreetmap.org/relation/3718820#map=7/49.472/-119.286
Should stretch all the way to Glacier NP I think
Looks to be one line of code with gdal. Then maybe try vectorizing the result? Then it can be put into mvt too
I forget where I am with US forest topo maps. But in any case I need to figure out the gdal commands to cut them into seamless tiles. Do I need to create a metadata.json file for these raster tiles?
Figure out why current tracks are shorter than the "known" length of the trail (2653mi or 2662mi to Manning Park). The distances I'm getting are 2D: 2592 mi, 3D: 2607.9 mi.
hm = Halfmile()
trail = hm.trail_full(alternates=False)
merged = linemerge([*trail.geometry])
projected = reproject(merged, geom.WGS84, geom.CA_ALBERS)
projected.length
# 4171048.0624210313 m
from haversine import haversine
dists = []
for a, b in zip(merged.coords, merged.coords[1:]):
dist = haversine((a[1], a[0]), (b[1], b[0]))
dists.append(dist)
sum(dists)
# 4171.3919978317 km
dists = []
for a, b in zip(projected.coords, projected.coords[1:]):
x = b[0] - a[0]
y = b[1] - a[1]
z = b[2] - a[2]
dist_m = sqrt(pow(x, 2) + pow(y, 2) + pow(z, 2))
dists.append(dist_m)
sum(dists)
# 4196977.125911263 m
Now that I've figured out how to use (online) vector tiles as an overlay/Fill layer in the mobile app (nst-guide/mobile#9), create vector tiles with (all?) overlays.
I think just an identifier is best because then you never have to update these vector tiles.
Because the individual layer file sizes are so small, you should probably package as many of these together as you can. So have basically three vector tilesets:
Because the file sizes are small, package up lots of layers into a single non-default tileset.
TODO
comment in 19d591c. It's been assigned to @kylebarron because they committed the code.I.e. use a DEM and find all the areas above 10,000 feet in Yosemite, or whatever elevation it is. Do that for each area
e.g. for national forests, national parks, wildernesses, also include a point layer for each that has its label at the polygon's pole of inaccessibility. You can use shapely to find the pole of inaccessibility
https://github.com/Toblerity/Shapely/blob/master/shapely/algorithms/polylabel.py
Note how to install conda and conda env create -f environment.yml
TODO
comment in 15eee87. It's been assigned to @kylebarron because they committed the code.TODO
comment in 15eee87. It's been assigned to @kylebarron because they committed the code.Now I've cleaned up the original transit api code. The question is how to format this to be consumed by the website.
There are naturally two layers, a route/line layer and a stop/point layer. The route should at a minimum have properties for
The stop should have at a minimum properties for:
Things to validate:
["get", "color"]
, and it's formatted correctly as rgbSee how difficult it would be to port the TNM style template to a Mapbox gl style:
https://viewer.nationalmap.gov/tools/topotemplate/
TODO
comment in 67ae82c. It's been assigned to @kylebarron because they committed the code.I.e. I have a National Parks boundaries dataset. What are interesting ways to connect that data to the trail?
Attributes:
When intersecting wilderness areas with the trail, wilderness areas may be a single MultiPolygon for a wilderness area. So when I intersect wilderness areas with the trail, I see polygons way out in Nevada, that aren't connected to California.
Instead, first split the multipolygons into individual polygons (keeping the attributes), then intersect. It might even be better to do two steps, first an intersection with the multipolygons, and then split into polygons, because then there's less data in the reshape.
data/code/data_source/transit_land.py
Lines 12 to 17 in 67f79b0
TODO
comment in 67f79b0. It's been assigned to @kylebarron because they committed the code.The missing block there is: 485212052_Skagit_Peak_FSTopo.tif
from this page: https://data.fs.usda.gov/geodata/rastergateway/states-regions/quad-index.php?blockID=48120
It's not created by get_quads
, so it's an error in my code
quad shapefiles: https://www.umesc.usgs.gov/data_library/maps_quads_figs/quad_indexes.html
Thoughts for generating raster tiles
I want to generate a list of tiles that's within x
miles from the trail. Using a GeoJSON file of a buffer polygon, use supermercado
to do:
cat buffer.geojson | supermercado burn 10
If you want to see the tiles on geojson.io, pass the list of tiles back to mercantile
:
cat buffer.geojson | supermercado burn 10 | mercantile shapes --collect > tiles.geojson
gdal2tiles.py
: https://gdal.org/programs/gdal2tiles.html
https://github.com/vincentsarago/lambda-tiler
Supermercado: https://github.com/mapbox/supermercado
Mercantile: https://github.com/mapbox/mercantile
Take two miles around trail, inside just Yosemite to start. find all water body features in that buffer, generate 100ft buffer.
Also 4 mile buffer around tuoulmne Meadows and Yosemite valley.
TODO
comment in e7b1afc. It's been assigned to @kylebarron because they committed the code.I tested out https://github.com/mapbox/vtshaver on the current OSM Liberty style.json, and it took 14/2525/*.pbf from 328KB to 139KB.
The size of the whole 14/**/*.pbf zoom level for CA, OR, and WA is 2.93GB, so if this ratio held for the entire area, that could go down to 1.2GB! It would be nothing to download everything within 20 miles of the trail.
TODO
comment in 164ca9c. It's been assigned to @kylebarron because they committed the code.Imagery is collected in the Military Grid Reference System. Here's a KML you can use to find imagery intersections with the trail.
The Sentinel AWS Readme lists how the file structure is set up. So for example, the UTM grid square that intersects Glacier Peak in WA, on July 22, 2019 is here:
> aws s3 ls s3://sentinel-s2-l2a/tiles/10/U/FU/2019/7/22/0/ --request-payer
PRE R10m/
PRE R20m/
PRE R60m/
PRE auxiliary/
PRE qi/
2019-07-22 21:15:40 553239 metadata.xml
2019-07-22 21:15:49 1035 productInfo.json
2019-07-22 21:15:40 1491 tileInfo.json
The tileInfo.json
has helpful information, including a GeoJSON representation of the geometry of the square, (not in WGS84), cloudyPixelPercentage
(!!), and the path to the image.
{
"path" : "tiles/10/U/FU/2019/7/22/0",
"timestamp" : "2019-07-22T19:11:13.752Z",
"utmZone" : 10,
"latitudeBand" : "U",
"gridSquare" : "FU",
"datastrip" : {
"id" : "S2A_OPER_MSI_L2A_DS_SGS__20190722T231540_S20190722T190551_N02.13",
"path" : "products/2019/7/22/S2A_MSIL2A_20190722T185921_N0213_R013_T10UFU_20190722T231540/datastrip/0"
},
"tileGeometry" : {
"type" : "Polygon",
"crs" : {
"type" : "name",
"properties" : {
"name" : "urn:ogc:def:crs:EPSG:8.8.1:32610"
}
},
"coordinates" : [ [ [ 600000.0, 5400000.0 ], [ 709800.0, 5400000.0 ], [ 709800.0, 5290200.0 ], [ 600000.0, 5290200.0 ], [ 600000.0, 5400000.0 ] ] ]
},
"tileDataGeometry" : {
"type" : "Polygon",
"crs" : {
"type" : "name",
"properties" : {
"name" : "urn:ogc:def:crs:EPSG:8.8.1:32610"
}
},
"coordinates" : [ [ [ 600001.0, 5399999.0 ], [ 709799.0, 5399999.0 ], [ 709799.0, 5290201.0 ], [ 600001.0, 5290201.0 ], [ 600001.0, 5399999.0 ] ] ]
},
"tileOrigin" : {
"type" : "Point",
"crs" : {
"type" : "name",
"properties" : {
"name" : "urn:ogc:def:crs:EPSG:8.8.1:32610"
}
},
"coordinates" : [ 600000.0, 5400000.0 ]
},
"dataCoveragePercentage" : 100.0,
"cloudyPixelPercentage" : 0.0,
"productName" : "S2A_MSIL2A_20190722T185921_N0213_R013_T10UFU_20190722T231540",
"productPath" : "products/2019/7/22/S2A_MSIL2A_20190722T185921_N0213_R013_T10UFU_20190722T231540"
}
If you look into a given folder for (utm square / year / month / day / snapshot number / resolution), you'll see this:
> aws s3 ls s3://sentinel-s2-l2a/tiles/10/U/EU/2019/7/22/0/R60m/ --request-payer
2019-07-22 21:07:22 512645 AOT.jp2
2019-07-22 21:07:22 2840382 B01.jp2
2019-07-22 21:07:22 3182126 B02.jp2
2019-07-22 21:07:22 3344222 B03.jp2
2019-07-22 21:07:22 3297405 B04.jp2
2019-07-22 21:07:22 3529444 B05.jp2
2019-07-22 21:07:22 3755698 B06.jp2
2019-07-22 21:07:22 3736390 B07.jp2
2019-07-22 21:07:22 3341917 B08.jp2
2019-07-22 21:07:22 3750817 B09.jp2
2019-07-22 21:07:22 3707894 B11.jp2
2019-07-22 21:07:22 3553773 B12.jp2
2019-07-22 21:07:22 3727852 B8A.jp2
2019-07-22 21:07:22 317802 SCL.jp2
2019-07-22 21:07:22 3711655 TCI.jp2
2019-07-22 21:07:22 3152856 WVP.jp2
Each of the B*.jp2 files are a different band. Wikipedia explains which band is which. However, easiest for my purposes is TCI.jp2
, which stands for True Color Image I think. For the 10m resolution, each file is around 129MB.
If you download a few images and run gdal2tiles.py --srcnodata="0,0,0" *.jp2 tiles
, you can get an image like this:
--srcnodata="0,0,0"
means that all fully black areas will become transparent, as seen above. The path of the satellites means that it won't always be simple to find continuous bands to display.
Note that Lambda isn't really a great place for this processing, since it has a scratch disk space limit of 500MB.
Note that while signed in to the aws cli, you need to add --request-payer
to do anything with a requester pays bucket. I.e.
aws s3 ls s3://sentinel-s2-l2a/ --request-payer
these values are in meters
total elevation change is: 346781.89000000095
total distance is: 4171048.0624210313
average slope is: 0.08314022874114661
from math import sqrt
from shapely.ops import linemerge
import geom
from data_source import Halfmile
def main():
hm = Halfmile()
gdf = hm.trail_full(False)
# Reproject to meters
gdf = geom.reproject_gdf(gdf, geom.WGS84, geom.CA_ALBERS)
line = linemerge(gdf.geometry.values)
# Calculate horizontal distance and vertical distance for each segment of
# the data:
dist_list, dz_list = calculate_dist_and_elev_change(line)
total_elevation_change = sum([abs(x) for x in dz_list])
total_distance = sum(dist_list)
average_slope = total_elevation_change / total_distance
print(f'total elevation change is: {total_elevation_change}')
print(f'total distance is: {total_distance}')
print(f'average slope is: {average_slope}')
def calculate_dist_and_elev_change(line):
dist_list = []
dz_list = []
for coord, next_coord in zip(line.coords, line.coords[1:]):
dx = next_coord[0] - coord[0]
dy = next_coord[1] - coord[1]
dz = next_coord[2] - coord[2]
dist = sqrt(dx ** 2 + dy ** 2)
dist_list.append(dist)
dz_list.append(dz)
return dist_list, dz_list
if __name__ == '__main__':
main()
TODO
comment in f26b9a6. It's been assigned to @kylebarron because they committed the code.TODO
comment in 3a64106. It's been assigned to @kylebarron because they committed the code.TODO
comment in dc2bd43. It's been assigned to @kylebarron because they committed the code.TODO
comment in 1e17d36. It's been assigned to @kylebarron because they committed the code.I.e. take all iPhone photos that have GPS information on that day, and try to use those to geotag that day's camera photos
Lines 57 to 60 in 8e37b6e
TODO
comment in 8e37b6e. It's been assigned to @kylebarron because they committed the code.TODO
comment in 0e11b44. It's been assigned to @kylebarron because they committed the code.TODO
comment in 660e995. It's been assigned to @kylebarron because they committed the code.Also check out labeling. Would be great to add labels like I currently have on streams to each section and each PCT alternate.
TODO
comment in 164ca9c. It's been assigned to @kylebarron because they committed the code.I need to figure out which x,y,z tile coordinates are within a polygon, I.e. within a given distance of the trail.
You can probably just use supermercado for this. https://github.com/mapbox/supermercado/blob/master/README.md
I think mercantile only finds bounding box intersections and not for complex polygons
This can definitely be a larger buffer than normal. Say for a buffer of like 10 miles around the trail, find every osm point (just point?) that's tagged as a peak or pass. Many of these will already be linked to Wikipedia, but for the others you can do a geosearch at that point, for that title, with a radius of like 100-300m.
For example this page doesn't exist:
https://inciweb.nwcg.gov/incident/4095
TODO
comment in d681daa. It's been assigned to @kylebarron because they committed the code.A declarative, efficient, and flexible JavaScript library for building user interfaces.
๐ Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. ๐๐๐
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google โค๏ธ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.