On 10/14/22 10:39, Stefan Gofferje wrote:
After the relative ease of my first mapserver project, I'm now motivated
and keen on getting another project running.
Since about a year (plus minnus), I'm every day checking for Sentinel I
and II data covering the area where I live. If there is anything, I
download the data and run some scripts which do automatic processing and
generating of composites, like SAR VV/VH or NDVI, NDCI, etc.
The ready images are then posted to a Telegram channel.
At the moment, there's still a number of manual steps (search data,
decide what to download, trigger download, trigger processing and
trigger Telegram push) and the existing scripts are Bash scripts calling
GDAL tools.
My first step now is to learn Python and move everything from Bash to
Python and automate it, ideally so far that the whole process can run in
a Docker container as a cron job.
Then I was thinking, instead of pushing images to Telegram, I could use
mapserver to serve the images to some web frontend. The map files are
very easily created programmatically, so that should be easy. I'm
missing a few things for that, though:
While going through mapserver's docs, I noticed that it can handle time
dimensions for imagery. That sounds like something worth exploring. The
docs describe the use of a shape file but unfortunately, I have no clue
how to to *create* that programmatically.
Regarding the frontend, I could just put up a single page leaflet app.
That would be easy enough. But I'm wondering if there's already some
cool solutions out there.
So, If anybody has a link or pointer for me on the programmatic creation
of shape files for time-indexed imagery and/or interesting web frontend
projects, I would be very grateful!
In case it is of use for anybody, here is the solution I came up with
last night:
1.) Add a time index metadatum to the output file
(gdaltranslate/gdal_edit -mo TIMEINDEX="-mm-ddThhmmZ"
2.) Hack together a Python script:
#!/usr/bin/python3
import os, sys
from osgeo import gdal
import geopandas as gpd
from shapely.geometry import box
import warnings
warnings.filterwarnings("ignore")
StartDir = str(sys.argv[1])
def getBounds(path):
raster = gdal.Open(path)
ulx, xres, xskew, uly, yskew, yres = raster.GetGeoTransform()
lrx = ulx + (raster.RasterXSize * xres)
lry = uly + (raster.RasterYSize * yres)
return box(lrx, lry, ulx, uly)
df = gpd.GeoDataFrame(columns=['location', 'geometry','timestamp'])
for dir, subdir, files in os.walk(StartDir):
for fname in files:
if fname.endswith(".tif"):
fullname = os.path.join(dir+"/", fname)
print (fullname)
ds=gdal.Open(fullname)
metadata=ds.GetMetadata()
ds=None
print(metadata)
df = df.append({'location': fname, 'geometry':
getBounds(fullname),'timestamp': metadata['TIMESTAMP']}, ignore_index=True)
#df = gpd.pd.concat(df,{'location': fname, 'geometry':
getBounds(os.path.join(dir+"/", fname))}, ignore_index=True)
df.to_file("tile-index.shp")
NOTE: suppressing all warnings is temporary until I figured out how to
move from df.append to pd.concat.
-Stefan
--
(o_ Stefan Gofferje| SCLT, MCP, CCSA
//\ Reg'd Linux User #247167 | VCP #2263
V_/_ https://www.gofferje.net | https://www.saakeskus.fi
___
MapServer-users mailing list
MapServer-users@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/mapserver-users
