Shakemap issueshttps://git.gfz-potsdam.de/groups/shakemap/-/issues2021-04-28T11:54:16+02:00https://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/2Add CI/CD pipeline2021-04-28T11:54:16+02:00Marius Kriegerowskimarius@gfz-potsdam.deAdd CI/CD pipelineWe just setup a CI/CD server for our group. We can recycle that until we have the need for more computational power.We just setup a CI/CD server for our group. We can recycle that until we have the need for more computational power.Marius Kriegerowskimarius@gfz-potsdam.deMarius Kriegerowskimarius@gfz-potsdam.dehttps://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/3Tools and standards2021-01-28T17:36:17+01:00Marius Kriegerowskimarius@gfz-potsdam.deTools and standardsI think it would be good to enforce a set of rules and tools right from the start to achieve a high code quality and standardization.
I'm open to other testing frameworks but I found `pytest` the best framework especially when it comes ...I think it would be good to enforce a set of rules and tools right from the start to achieve a high code quality and standardization.
I'm open to other testing frameworks but I found `pytest` the best framework especially when it comes to testing APIs as it comes with a great set of features and is simpler than using python's unittest framework. I will disable pushing to master branch so we only merge into that after approval of at least one reviewer.
Since recently I wasn't a fan of automatic code formatting myself but I was taught better. It helps tremendously to have common styling standards. Thus I would like to enforce style checks also in CI using `pylint` and `black`.
Opinions? Other tools you would like to start with right from the beginning @ds @gweather @rizac?https://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/5add top level license file2021-02-16T10:19:08+01:00Marius Kriegerowskimarius@gfz-potsdam.deadd top level license filehttps://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/6extend .gitignore2021-02-16T10:18:46+01:00Marius Kriegerowskimarius@gfz-potsdam.deextend .gitignoreMaybe use auto generated?Maybe use auto generated?https://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/7Proposal simple Regionalization storage2021-04-09T14:55:12+02:00Riccardo ZaccarelliProposal simple Regionalization storageFollowing our discussion, here a sketch of what is roughly implemented in eGSIM. The main idea is to try the same here and see how to proceed afterwards, e.g., move to geopandas, or remove the eGSIM implementation and make this an import...Following our discussion, here a sketch of what is roughly implemented in eGSIM. The main idea is to try the same here and see how to proceed afterwards, e.g., move to geopandas, or remove the eGSIM implementation and make this an imported library of eGSIM, although this is maybe not really feasible because eGSIM works on a SQL database behind the hood (sorry I did not think about it in the call).
Anyway, a simple regionalization can be achieved via two different JSON files, a regionalization and a "gsim selection".
1. The regionalization is a [GeoJSON](https://en.wikipedia.org/wiki/GeoJSON) object, basically a list of geographic regions, with associated tectonic region type (trt):
```json
{
"type": "FeatureCollection",
"features": [
{
"type": "Feature",
"geometry": {
"type": "Polygon",
"coordinates": [[102.0, 0.5], ... ]
},
"properties": {
"trt": "Active shallow crust"
}
},
...
]
}
```
2. A a "gsim selection" is a JSON object, i.e. a `dict` mapping a trt to a list of GSIMs (just an example with the first two random GSIMs I can recall):
```json
{
"Active shallow crust": ["AkkarBommer2010", "CauzziFaccioli2008"],
...
}
```
In eGSIM, given a point, we iterate over each region in the geoJSON and check if the point is in the region (using the methods of the Python `shapely` library [I guess](https://github.com/rizac/eGSIM/blob/2836bdfe6af84cf9e5d96b97f401390576d96cbc/egsim/models.py#L371)). When a point is in the region, we return the GSIMs of that region from the "gsim selection" JSON object.
From my experience, the only issue I had with this approach is naming conventions (different TRTs denoting the same concept), which is solved in eGSIM with another database table (here we might define our unique names and be strict). Obviously, there is also scalability: for hundreds of regions I never had problems. Above, I don't know.
Last note: in eGSIM we have several "gsim selection" dicts, one for each research project model, e.g. SHARE, ESHM2020 (still ongoing), but this is not an issue here I guess, so we can stick to a single "gsim selection" dicthttps://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/8Migrate "Event" Object from Shakeyground (rename to Earthquake)2021-02-25T17:04:26+01:00Graeme WeatherillMigrate "Event" Object from Shakeyground (rename to Earthquake)Creation of the initial earthquake object: refactor the previous definition from the shakyground code, extend methods for rupture generation, add tests and documentationCreation of the initial earthquake object: refactor the previous definition from the shakyground code, extend methods for rupture generation, add tests and documentationGraeme WeatherillGraeme Weatherillhttps://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/9Migrate "Synthetic Rupture Generator from Shakyground2021-03-11T16:43:34+01:00Graeme WeatherillMigrate "Synthetic Rupture Generator from ShakygroundMigrate the synthetic rupture generator from the previous shakyground code. Most of the original code can be ported without change, but test coverage should be improved and documentation added.Migrate the synthetic rupture generator from the previous shakyground code. Most of the original code can be ported without change, but test coverage should be improved and documentation added.Graeme WeatherillGraeme Weatherillhttps://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/10Migrate simple Shakemap object from shakyground2021-03-11T16:43:22+01:00Graeme WeatherillMigrate simple Shakemap object from shakygroundMigrate and the basic object for generating and storing the shakemaps and refactor using recent OpenQuake changes. Takes as input an event object, a set of sites and a dictionary of GMPEs, and returns the ground motion values for the sit...Migrate and the basic object for generating and storing the shakemaps and refactor using recent OpenQuake changes. Takes as input an event object, a set of sites and a dictionary of GMPEs, and returns the ground motion values for the sites. Ground motion fields for the median ("average") field and the standard deviations. Requires tests and improved documentation.Graeme WeatherillGraeme Weatherillhttps://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/12Add Site Model Object2021-02-25T17:04:34+01:00Graeme WeatherillAdd Site Model ObjectThe handling of the target sites for calculating the ground motions can be improved with respect to the previous implementation. Create a general site model object to control the I/O of the target site information and scientific function...The handling of the target sites for calculating the ground motions can be improved with respect to the previous implementation. Create a general site model object to control the I/O of the target site information and scientific functionalities needed to characterise the sites for the selected ground motion model. Explore inheriting/extending OpenQuake SiteCollection object, otherwise via composition.Graeme WeatherillGraeme Weatherillhttps://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/15Accept/parse QuakeML from GEOFON to produce input for a new run2021-05-20T14:45:48+02:00Peter EvansAccept/parse QuakeML from GEOFON to produce input for a new runWhen a new event occurs, we can query the [GEOFON](http://geofon.gfz-potsdam.de/eqinfo/list.php) [fdsnws-event web service](http://geofon.gfz-potsdam.de/fdsnws/event/1) to get event parameters, including focal mechanism if available.
An...When a new event occurs, we can query the [GEOFON](http://geofon.gfz-potsdam.de/eqinfo/list.php) [fdsnws-event web service](http://geofon.gfz-potsdam.de/fdsnws/event/1) to get event parameters, including focal mechanism if available.
Angelo:
> Peter to work with Riccardo in preparing the input which will be as simple as possible for
the moment (yaml/json) including EQ Lat/Lon/Depth/Magnitude, and as soon as available
the focal mechanism
This merge request does the query, parses the XML, and produces a basic Python structure. This can be easily serialised in the flavour of your choice.
For example, here is the YAML output for a large event in Greece: [gfz2021eksc](http://geofon.gfz-potsdam.de/eqinfo/event.php?id=gfz2021eksc):
```
date: '2021-03-04'
focalmechanism:
- dip: 48.75643809
rake: -91.74289876
strike: 146.5577511
- dip: 41.27378495
rake: -88.01325753
strike: 329.2004058
id: gfz2021eksc
magnitude: 6.278409404
origin:
depth: 10.0
latitude: 39.79063416
longitude: 22.28781891
preferred_focalmechanism_id: smi:org.gfz-potsdam.de/geofon/FocalMechanism/20210304185910.015761.20
39019
preferred_magnitude_id: smi:org.gfz-potsdam.de/geofon/Magnitude/20210304185910.01602.2039021
preferred_origin_id: smi:org.gfz-potsdam.de/geofon/Origin/20210304192045.9955.2520068
time: 18:38:20.19439Z
```
There are two (dip, rake, strike) sections for the two nodal planes. The QuakeML doesn't help us decide which to use.
From here, some integration will be needed to use the YAML to drive a run of the shake map.https://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/16Add GeoTiff and ARC/INFO ASCII GRID exporter for shakemap (as raster)2021-04-28T11:04:45+02:00Graeme WeatherillAdd GeoTiff and ARC/INFO ASCII GRID exporter for shakemap (as raster)For integration in the EarthExplorer platform (amongst other users) exporters are needed to write the shakemaps into common geospatial data formats for visualisation. Two formats should be considered initially:
1. GeoTiff - The common a...For integration in the EarthExplorer platform (amongst other users) exporters are needed to write the shakemaps into common geospatial data formats for visualisation. Two formats should be considered initially:
1. GeoTiff - The common and compact standard for single band raster data
2. ARC/INFO ASCII GRID - Relatively simple text format recognised by all major GIS platforms. Easy to build using numpy alone, though less compact than GeoTiff.
For handling GeoTiff formats I usually use `rasterio` (https://rasterio.readthedocs.io/en/latest/), which brings few dependencies and is easy to install and run. It can read and write most common raster formats. Another option is to use the Python GDAL bindings (https://pypi.org/project/GDAL/). This is a more powerful tool and would potentially open up the possibility to write to a wider variety of raster formats in the future. However, it creates a dependency on GDAL, which in my experience has been extremely challenging to run on a lot of platforms (I have only ever had bad experiences getting it to work cleanly, if at all, on OSX, and run into problems on both Windows and Linux).
@marius @rizac @nils @eggi Any suggestions on the most suitable Python tools for building rasters?https://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/18Shakemaps failing when crossing the international date line (IDL)2021-06-18T10:47:46+02:00Graeme WeatherillShakemaps failing when crossing the international date line (IDL)When an event occurs close to the international dateline (180˚E/˚W) the target site model and resulting shakemap cannot be built.
There are several steps to solving this problem:
1. OpenQuake site models require longitudes in the range...When an event occurs close to the international dateline (180˚E/˚W) the target site model and resulting shakemap cannot be built.
There are several steps to solving this problem:
1. OpenQuake site models require longitudes in the range -180˚E to 180˚E, so permit bounding boxes to generate limits beyond -180˚E and 180˚E in order to generate the site models, but then translate sites in the range long < -180˚E and > 180˚E accordingly. This would make it possible to calculate the ground motions for the sites.
2. For exporting to raster there are two options: i) Split the raster into two different rasters (one either side of the IDL), ii) translate the sites back so that longitudes can be outside of the -180˚E to 180˚E band and retain the simple bounding box shape.
@nils and @eggi - what is the preferred approach for EarthExplorer in this case? Split the raster into two rasters or tolerate longitude values < -180˚E and > 180˚Ehttps://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/20Add shakemap summary statistics into ShakemapsWorkflowResult2021-06-02T14:23:01+02:00Graeme WeatherillAdd shakemap summary statistics into ShakemapsWorkflowResultFor the case where a user may be building an archive of shakemaps (e.g. EarthquakeExplorer) it can be useful to have some summary statistics about the shakemaps (e.g. mean acceleration, max acceleration etc.) included into the `Shakemaps...For the case where a user may be building an archive of shakemaps (e.g. EarthquakeExplorer) it can be useful to have some summary statistics about the shakemaps (e.g. mean acceleration, max acceleration etc.) included into the `ShakemapsWorkflowResult` object. This could make it easier to store attributes that could be used for filtering/searching the shakemaps subsequently, rather than having to retrieve the values from the rastershttps://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/21Add shakemaps for macroseismic intensity2021-06-15T08:55:48+02:00Graeme WeatherillAdd shakemaps for macroseismic intensityNeeded for i) EarthquakeExplorer, ii) use cases where the shakemap may be linked to building damage or vulnerability calculations based on macroseismic intensity.
There are two options for doing this:
1) Define a global regionalization...Needed for i) EarthquakeExplorer, ii) use cases where the shakemap may be linked to building damage or vulnerability calculations based on macroseismic intensity.
There are two options for doing this:
1) Define a global regionalization for intensity prediction equations (IPEs)
2) Use a ground motion intensity conversion equation (GMICE)
Option 1 would be preferable if there were a sufficiently large number of intensity prediction equations to cover all tectonic region types in order to define these everywhere. In reality, there are considerably fewer IPEs in the literature than GMPEs, some are very localised (to a specific country or region) and not all cover all of the different tectonic environments needed. Therefore the alternative is to use the GMICEs. There are several available in the literature (e.g. Worden et al., 2012; Faenza & Michelini, 2010; Caprio et al., 2015) and though there are some regional differences in GMICEs, these are more appropriate for application worldwide.
GMICEs usually describe the relationship between PGA, PGV and/or SA and macroseismic intensity. Though PGV is often seen as the better predictor of intensity, not all GMPEs selected in the default regionalisation actually define PGV, therefore it would need to be inferred from PGA or SA. Given that this conversion would carry larger uncertainty than using PGA - intensity, I would stick with just the PGA to intensity conversions.https://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/22Datetime parsing in shakemaps_from_quakeml2021-05-27T12:55:10+02:00Nils BrinckmannDatetime parsing in shakemaps_from_quakemlCurrently there is a problem on processing some of the events from geofon with the `shakemaps_from_quakeml` workflow.
This effects for example the following event:
- gfz2021kdvq
- gfz2021keus
Traceback is the following:
```
Traceback ...Currently there is a problem on processing some of the events from geofon with the `shakemaps_from_quakeml` workflow.
This effects for example the following event:
- gfz2021kdvq
- gfz2021keus
Traceback is the following:
```
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/usr/src/app/shakyground2/shakyground2/workflows.py", line 93, in shakemaps_from_quakeml
results.earthquake = Earthquake.from_quakeml(event_id)
File "/usr/src/app/shakyground2/shakyground2/earthquake.py", line 219, in from_quakeml
d_t = datetime.datetime.fromisoformat(" ".join([event["date"], event_time]))
ValueError: Invalid isoformat string: '2021-05-25 09:41:17.20'
```https://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/23shakemaps_from_quakeml error on fetching event2021-05-27T12:55:57+02:00Nils Brinckmannshakemaps_from_quakeml error on fetching eventThere is an xml error:
```
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/usr/src/app/shakyground2/shakyground2/workflows.py", line 93, in shakemaps_from_quakeml
results.earthquake = Earthquake.fro...There is an xml error:
```
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/usr/src/app/shakyground2/shakyground2/workflows.py", line 93, in shakemaps_from_quakeml
results.earthquake = Earthquake.from_quakeml(event_id)
File "/usr/src/app/shakyground2/shakyground2/earthquake.py", line 209, in from_quakeml
event = fetch_quakeml(path)
File "/usr/src/app/shakyground2/shakyground2/io/import_fdsnws_eq.py", line 182, in fetch_quakeml
root = ET.fromstring(buf)
File "/usr/local/lib/python3.8/xml/etree/ElementTree.py", line 1321, in XML
return parser.close()
xml.etree.ElementTree.ParseError: no element found: line 1, column 0
```
for the following event:
- gfz2021jwsm
Code to reproduce:
```
from shakyground2.workflows import shakemaps_from_quakeml
results = shakemaps_from_quakeml("gfz2021jwsm")
```https://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/26Tiff output coordinate system problem2021-12-02T12:54:04+01:00Nils BrinckmannTiff output coordinate system problemWhen I run `gdalinfo` over a tiff that comes out of the `shakemaps_from_quakeml` workflow the coordinate system as a little problem:
```
Driver: GTiff/GeoTIFF
Files: PGA.tiff
Size is 2107, 1740
Coordinate System is:
GEOGCS["WGS 84",
...When I run `gdalinfo` over a tiff that comes out of the `shakemaps_from_quakeml` workflow the coordinate system as a little problem:
```
Driver: GTiff/GeoTIFF
Files: PGA.tiff
Size is 2107, 1740
Coordinate System is:
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.257223563,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]]
Origin = (89.594909899532766,27.287898476052039)
Pixel Size = (0.008333333333333,0.008333333333333)
Metadata:
AREA_OR_POINT=Area
Image Structure Metadata:
COMPRESSION=DEFLATE
INTERLEAVE=BAND
Corner Coordinates:
Upper Left ( 89.5949099, 27.2878985) ( 89d35'41.68"E, 27d17'16.43"N)
Lower Left ( 89.5949099, 41.7878985) ( 89d35'41.68"E, 41d47'16.43"N)
Upper Right ( 107.1532432, 27.2878985) (107d 9'11.68"E, 27d17'16.43"N)
Lower Right ( 107.1532432, 41.7878985) (107d 9'11.68"E, 41d47'16.43"N)
Center ( 98.3740766, 34.5378985) ( 98d22'26.68"E, 34d32'16.43"N)
Band 1 Block=2107x1 Type=Float64, ColorInterp=Gray
```
The corner coordinates upper and lower values are reversed.
It is possible to fix those with
```
gdalwarp PGA.tiff PGA_fixed.tiff -of GTiff -s_srs '+proj=latlong' -t_srs 'epsg:4326' -overwrite
```
but it would be way better if we could give back a tiff with the right corner coordinates.
The actual spatial projection of the single pixels to points on the globe is not effected by this issue - it seems that just the ordering of the values is different then expected by gdal.
One of the reasons to fix it, is to allow the usage of the tiffs in gmt grdimage. Non fixed versions can result in an error when allocating memory.https://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/27Fix IQR statistic for shakemaps_from_quakeml workflow2022-01-19T16:07:06+01:00Nils BrinckmannFix IQR statistic for shakemaps_from_quakeml workflowWhile the min, max & median values in the statistics dict are as fractiles of `g`, the IQR value is not.
See: https://git.gfz-potsdam.de/shakemap/shakyground2/-/blob/master/shakyground2/workflows.py#L261
I would expect something like:
...While the min, max & median values in the statistics dict are as fractiles of `g`, the IQR value is not.
See: https://git.gfz-potsdam.de/shakemap/shakyground2/-/blob/master/shakyground2/workflows.py#L261
I would expect something like:
```
"IQR": np.exp(np.quantile(mean_shakemaps[imt], 0.75)) - np.exp(np.quantile(mean_shakemaps[imt], 0.25))
```https://git.gfz-potsdam.de/shakemap/shakyground2/-/issues/28Bug in shakemaps_from_quakeml workflow for gfz2022gcsb2022-03-31T11:28:00+02:00Nils BrinckmannBug in shakemaps_from_quakeml workflow for gfz2022gcsbWhen I run the following code:
```python
from shakyground2.workflows import shakemaps_from_quakeml
imts = ["PGA", "SA(0.3)", "SA(1.0)", "MMI"]
event_id = "gfz2022gcsb"
workflow_results = shakemaps_from_quakeml(event_id=event_id, imts=...When I run the following code:
```python
from shakyground2.workflows import shakemaps_from_quakeml
imts = ["PGA", "SA(0.3)", "SA(1.0)", "MMI"]
event_id = "gfz2022gcsb"
workflow_results = shakemaps_from_quakeml(event_id=event_id, imts=imts, config=None)
```
Then I get an error:
```
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/usr/src/app/shakyground2/shakyground2/workflows.py", line 232, in shakemaps_from_quakeml
results.contours["mean"][imt] = shakemap.get_contours(
File "/usr/src/app/shakyground2/shakyground2/shakemap.py", line 542, in get_contours
geometries.append(LineString([(row[0], row[1]) for row in seg]))
File "/usr/local/lib/python3.8/dist-packages/shapely/geometry/linestring.py", line 46, in __init__
ret = geos_linestring_from_py(coordinates)
File "shapely/speedups/_speedups.pyx", line 156, in shapely.speedups._speedups.geos_linestring_from_py
ValueError: LineStrings must have at least 2 coordinate tuples
```
From my understanding it happens when the contour lines should be generated. Maybe it is assosiated with the used libraries that we use in the Eq Explorer.
See https://git.gfz-potsdam.de/earthquake-explorer/gfz-earthquake-explorer/-/issues/254 as well.