diff --git a/.github/workflows/engine_pr_test.yml b/.github/workflows/engine_pr_test.yml index 2ef3524cdc8e..cf3c562a9d80 100644 --- a/.github/workflows/engine_pr_test.yml +++ b/.github/workflows/engine_pr_test.yml @@ -70,6 +70,7 @@ jobs: ruff check --preview . pytest calculators -k 'risk or damage or bcr' -x --doctest-modules --disable-warnings --color=yes --durations=5 pytest -xs --doctest-modules --disable-warnings --color=yes --durations=8 sep hmtk risklib commonlib baselib hazardlib + pytest --doctest-modules ../doc/contributing/*.rst server_public_mode: runs-on: ${{ matrix.os }} diff --git a/debian/changelog b/debian/changelog index e3fb052b44b1..00499633430f 100644 --- a/debian/changelog +++ b/debian/changelog @@ -1,4 +1,6 @@ [Michele Simionato] + * Internal: changed the ordering in the composite risk model from + (loss_type, riskid) -> (riskid, loss_type) * Added an exporter for trt_gsim * Internal: added utility function readinput.read_source_models diff --git a/doc/contributing/developing-with-the-engine.rst b/doc/contributing/developing-with-the-engine.rst index 9bd94100e39c..597374567297 100644 --- a/doc/contributing/developing-with-the-engine.rst +++ b/doc/contributing/developing-with-the-engine.rst @@ -262,7 +262,8 @@ Then the hazard curve can be computed as follows:: >>> sitecol = readinput.get_site_collection(oq) >>> gsims = readinput.get_gsim_lt(oq).values['*'] >>> calc_hazard_curve(sitecol, src, gsims, oq) - [[0.00507997]]> + array([[0.00508004]], dtype=float32) + Working with GMPEs directly: the ContextMaker --------------------------------------------- @@ -345,6 +346,61 @@ the result provided by ``calc_hazard_curve(sitecol, src, gsims, oq)`` in the sec If you want to know exactly how ``get_pmap`` works you are invited to look at the source code in ``openquake.hazardlib.contexts``. +Generating ground motion fields from a rupture +---------------------------------------------- + +The easiest way to create a finite size rupture (a.k.a. planar rupture) +is to use the factory function `get_planar`: + +>>> from openquake.hazardlib.source.rupture import get_planar + +The function requires in input a site and a magnitude scaling relationship, +so first you have to build such objects: + +>>> [site] = sitecol # since there is a single site +>>> from openquake.hazardlib.scalerel import WC1994 +>>> msr = WC1994() # magnitude scaling relationship +>>> mag = 6. +>>> rup = get_planar(site, msr, mag, aratio=1., strike=11., dip=38., +... rake=55., trt=cmaker.trt) + +If you want to generate the GMF produced by a rupture (i.e. to emulate +a scenario calculation) you need to supplement the number of +occurrences of the rupture and a random seed, i.e. you need to convert the +hazardlib rupture into an EBRupture: + +>>> from openquake.hazardlib.source.rupture import EBRupture +>>> ebr = EBRupture(rup, n_occ=2, seed=42) + +Then you can use the GmfComputer class to perform the calculation: + +>>> from openquake.hazardlib.calc.gmf import GmfComputer +>>> gc = GmfComputer(ebr, sitecol, cmaker) +>>> gc.compute_all() # returns a DataFrame + gmv_0 eid sid rlz +0 0.660239 0 0 0 +1 0.301583 1 0 0 + +`gmv_0` is the value of the ground motion field for the first IMT (i.e PGA in this +case), `eid` the event ID, `sid` the site ID (there is a single site in this case) +and `rlz` the realization index. +In scenario calculations there is a realization for each GSIM and in this case +there is a single GSIM, so rlz=0. The total number of events is the number +of realizations times the number of occurrences and therefore in this case +the event ID (`eid`) can only have the values 0 or 1. + +It is also possible to perform calculations with point-like ruptures +(i.e. ignoring the finite-size effects): + +>>> from openquake.hazardlib.source.rupture import PointRupture +>>> occ_rate = None # not used in the GmfComputer +>>> rup = PointRupture(mag, cmaker.trt, site.location, occ_rate, cmaker.tom) +>>> ebr = EBRupture(rup, n_occ=2, seed=42) +>>> GmfComputer(ebr, sitecol, cmaker).compute_all() + gmv_0 eid sid rlz +0 0.541180 0 0 0 +1 0.247199 1 0 0 + Working with verification tables -------------------------------- @@ -399,10 +455,8 @@ Running the engine tests ------------------------ If you are a hazard scientist contributing a bug fix to a GMPE (or any other kind of bug fix) you may need to run the -engine tests and possibly change the expected files if there is a change in the numbers. The way to do it is to start -the dbserver and then run the tests from the repository root:: +engine tests and possibly change the expected files if there is a change in the numbers. The way to do it is to give the following command from the repository root:: - $ oq dbserver start $ pytest -vx openquake/calculators If you get an error like the following:: diff --git a/docker/scripts/oq-start.sh b/docker/scripts/oq-start.sh index be38e808d23b..4e60ec36dbc2 100755 --- a/docker/scripts/oq-start.sh +++ b/docker/scripts/oq-start.sh @@ -19,7 +19,7 @@ # This is required to load a custom local_settings.py when 'oq webui' is used. export PYTHONPATH=$HOME -oq engine --upgrade-db & +oq dbserver start & # Wait the DbServer to come up; may be replaced with a "oq dbserver wait" while : diff --git a/openquake/baselib/general.py b/openquake/baselib/general.py index 1d7b041da4fc..b43b61815a63 100644 --- a/openquake/baselib/general.py +++ b/openquake/baselib/general.py @@ -45,7 +45,7 @@ import numpy import pandas from decorator import decorator -from openquake.baselib import __version__ +from openquake.baselib import __version__, config from openquake.baselib.python3compat import decode U8 = numpy.uint8 @@ -436,7 +436,8 @@ def gettemp(content=None, dir=None, prefix="tmp", suffix="tmp", remove=True): if dir is not None: if not os.path.exists(dir): os.makedirs(dir) - fh, path = tempfile.mkstemp(dir=dir, prefix=prefix, suffix=suffix) + fh, path = tempfile.mkstemp(dir=dir or config.directory.custom_tmp, + prefix=prefix, suffix=suffix) if remove: _tmp_paths.append(path) with os.fdopen(fh, "wb") as fh: diff --git a/openquake/calculators/multi_risk.py b/openquake/calculators/multi_risk.py index 83274cc72439..b6cf65e721b6 100644 --- a/openquake/calculators/multi_risk.py +++ b/openquake/calculators/multi_risk.py @@ -58,7 +58,7 @@ def get_dmg_csq(crm, assets_by_site, gmf, time_event): for k, w in zip(df.risk_id, df.weight)] # NB: risk logic trees are not yet supported in multi_risk [peril] = rm.imt_by_lt.values() - fracs[li] = rm.scenario_damage(loss_type, assets, + fracs[li] = rm.scenario_damage('earthquake', loss_type, assets, pandas.DataFrame({peril: [gmv]})) csq = crm.compute_csq(assets, fracs, df, crm.oqparam) number = assets['value-number'] diff --git a/openquake/commands/compare.py b/openquake/commands/compare.py index 1aca44a8ab89..ff7e7ff3eb51 100644 --- a/openquake/commands/compare.py +++ b/openquake/commands/compare.py @@ -469,8 +469,8 @@ def compare_oqparam(calc_ids: int): ds1 = datastore.read(calc_ids[1]) dic0 = vars(ds0['oqparam']) dic1 = vars(ds1['oqparam']) - common = sorted(set(dic0) & set(dic1)) - for key in common: + common = set(dic0) & set(dic1) - {'hdf5path'} + for key in sorted(common): if dic0[key] != dic1[key]: print('%s: %s != %s' % (key, dic0[key], dic1[key])) diff --git a/openquake/commonlib/oqvalidation.py b/openquake/commonlib/oqvalidation.py index 594efe62e3b4..25ec02d52fc6 100644 --- a/openquake/commonlib/oqvalidation.py +++ b/openquake/commonlib/oqvalidation.py @@ -1472,12 +1472,12 @@ def set_loss_types(self): else: self._parent = None # set all_cost_types - # rt has the form 'vulnerability/structural', 'fragility/...', ... - costtypes = set(rt.rsplit('/')[1] for rt in self.risk_files) + # rt has the form 'earthquake/vulnerability/structural', ... + costtypes = set(rt.split('/')[2] for rt in self.risk_files) if not costtypes and self.hazard_calculation_id: try: self._risk_files = rfs = get_risk_files(self._parent.inputs) - costtypes = set(rt.rsplit('/')[1] for rt in rfs) + costtypes = set(rt.split('/')[2] for rt in rfs) except OSError: # FileNotFound for wrong hazard_calculation_id pass self.all_cost_types = sorted(costtypes) # including occupants diff --git a/openquake/engine/tests/aristotle1/aggrisk.txt b/openquake/engine/tests/aristotle1/aggrisk.txt index b0209a60ab59..e9c47bb9c40b 100644 --- a/openquake/engine/tests/aristotle1/aggrisk.txt +++ b/openquake/engine/tests/aristotle1/aggrisk.txt @@ -1,6 +1,6 @@ | gsim | weight | contents | nonstructural | structural | occupants | area | number | residents | |-------------------+--------+-------------+---------------+------------+-----------+--------+----------+-----------| -| [ChiouYoungs2014] | 0.3500 | 47_928_608 | 86_326_184 | 49_106_732 | 1.3142 | 31_790 | 107.8218 | 1_063 | -| [BooreEtAl2014] | 0.3500 | 71_780_840 | 76_948_504 | 55_768_280 | 0.4272 | 43_355 | 250.2 | 1_185 | -| [ZhaoEtAl2016Asc] | 0.3000 | 123_750_848 | 88_760_720 | 40_181_644 | 0.2262 | 61_196 | 349.9 | 1_304 | -| Average | 1.0000 | 79_023_561 | 83_774_357 | 48_760_747 | 0.6773 | 44_660 | 230.3 | 1_178 | \ No newline at end of file +| [ChiouYoungs2014] | 0.3500 | 68_505_824 | 77_741_872 | 48_760_420 | 1.3142 | 31_790 | 107.8218 | 1_063 | +| [BooreEtAl2014] | 0.3500 | 81_108_936 | 104_466_288 | 45_423_496 | 0.4272 | 43_355 | 250.2 | 1_185 | +| [ZhaoEtAl2016Asc] | 0.3000 | 105_540_112 | 68_223_000 | 39_220_012 | 0.2262 | 61_196 | 349.9 | 1_304 | +| Average | 1.0000 | 84_027_200 | 84_239_756 | 44_730_374 | 0.6773 | 44_660 | 230.3 | 1_178 | \ No newline at end of file diff --git a/openquake/hazardlib/contexts.py b/openquake/hazardlib/contexts.py index a5b66f36a2cf..7f95f594f774 100644 --- a/openquake/hazardlib/contexts.py +++ b/openquake/hazardlib/contexts.py @@ -199,13 +199,18 @@ class Oq(object): mea_tau_phi = False split_sources = True use_rates = False + af = None def __init__(self, **hparams): vars(self).update(hparams) @property def min_iml(self): - return numpy.array([1E-10 for imt in self.imtls]) + try: + imtls = self.imtls + except AttributeError: + imtls = self.hazard_imtls + return numpy.array([1E-10 for imt in imtls]) def get_reqv(self): if 'reqv' not in self.inputs: @@ -531,7 +536,10 @@ class ContextMaker(object): def __init__(self, trt, gsims, oq, monitor=Monitor(), extraparams=()): self.trt = trt - self.gsims = gsims + if isinstance(gsims, dict): + self.gsims = gsims + else: + self.gsims = {gsim: U32([i]) for i, gsim in enumerate(gsims)} # NB: the gid array can be overridden later on self.gid = numpy.arange(len(gsims), dtype=numpy.uint16) if isinstance(oq, dict): @@ -595,8 +603,6 @@ def _init1(self, param): self.disagg_bin_edges = param.get('disagg_bin_edges', {}) self.ps_grid_spacing = param.get('ps_grid_spacing') self.split_sources = self.oq.split_sources - - def _init2(self, param, extraparams): for gsim in self.gsims: if hasattr(gsim, 'set_tables'): if len(self.mags) == 0 and not is_modifiable(gsim): @@ -604,7 +610,8 @@ def _init2(self, param, extraparams): 'You must supply a list of magnitudes as 2-digit ' 'strings, like mags=["6.00", "6.10", "6.20"]') gsim.set_tables(self.mags, self.imtls) - self.effect = param.get('effect') + + def _init2(self, param, extraparams): for req in self.REQUIRES: reqset = set() for gsim in self.gsims: @@ -1221,7 +1228,9 @@ def get_4MN(self, ctxs, gsim): def get_att_curves(self, site, msr, mag, aratio=1., strike=0., dip=45., rake=-90): """ - :returns: 4 attenuation curves mu, sig, tau, phi + :returns: + 4 attenuation curves mea, sig, tau, phi + (up to 500 km from the site at steps of 5 km) """ from openquake.hazardlib.source import rupture rup = rupture.get_planar( diff --git a/openquake/hazardlib/source/point.py b/openquake/hazardlib/source/point.py index ca11408d41c8..e05c136e0e3f 100644 --- a/openquake/hazardlib/source/point.py +++ b/openquake/hazardlib/source/point.py @@ -305,8 +305,8 @@ def _gen_ruptures(self, shift_hypo=False, step=1, iruptures=False): for d, (drate, cdep) in hdd_[::step]: rate = mrate * nrate * drate yield PointRupture( - mag, np.rake, self.tectonic_region_type, - Point(clon, clat, cdep), np.strike, np.dip, rate, + mag, self.tectonic_region_type, + Point(clon, clat, cdep), rate, self.temporal_occurrence_model, self.lower_seismogenic_depth) diff --git a/openquake/hazardlib/source/rupture.py b/openquake/hazardlib/source/rupture.py index de156c6a6d2e..c5a02b771203 100644 --- a/openquake/hazardlib/source/rupture.py +++ b/openquake/hazardlib/source/rupture.py @@ -583,19 +583,15 @@ def get_cdppvalue(self, target, buf=1.0, delta=0.01, space=2.): class PointSurface: """ A fake surface used in PointRuptures. - The parameters `hypocenter` and `strike` are determined by - collapsing the corresponding parameters in the original PointSource. """ - def __init__(self, hypocenter, strike, dip): + def __init__(self, hypocenter): self.hypocenter = hypocenter - self.strike = strike - self.dip = dip def get_strike(self): - return self.strike + return 0 def get_dip(self): - return self.dip + return 0 def get_top_edge_depth(self): return self.hypocenter.depth @@ -628,17 +624,17 @@ class PointRupture(ParametricProbabilisticRupture): A rupture coming from a far away PointSource, so that the finite size effects can be neglected. """ - def __init__(self, mag, rake, tectonic_region_type, hypocenter, strike, - dip, occurrence_rate, temporal_occurrence_model, zbot): + def __init__(self, mag, tectonic_region_type, hypocenter, + occurrence_rate, temporal_occurrence_model, zbot=0): + self.mag = mag self.tectonic_region_type = tectonic_region_type self.hypocenter = hypocenter - self.mag = mag - self.strike = strike - self.rake = rake - self.dip = dip self.occurrence_rate = occurrence_rate self.temporal_occurrence_model = temporal_occurrence_model - self.surface = PointSurface(hypocenter, strike, dip) + self.surface = PointSurface(hypocenter) + self.rake = 0 + self.dip = 0 + self.strike = 0 self.zbot = zbot # bottom edge depth, used in Campbell-Bozorgnia @@ -719,7 +715,7 @@ class ExportedRupture(object): Simplified Rupture class with attributes rupid, events_by_ses, indices and others, used in export. - :param rupid: rupture.seed ID + :param rupid: rupture ID :param events_by_ses: dictionary ses_idx -> event records :param indices: site indices """ @@ -744,13 +740,14 @@ class EBRupture(object): """ seed = 'NA' # set by the engine - def __init__(self, rupture, source_id, trt_smr, n_occ=1, id=None, e0=0): + def __init__(self, rupture, source_id=0, trt_smr=0, n_occ=1, id=0, e0=0, seed=42): self.rupture = rupture self.source_id = source_id self.trt_smr = trt_smr self.n_occ = n_occ self.id = source_id * TWO30 + id self.e0 = e0 + self.seed = seed @property def tectonic_region_type(self): diff --git a/openquake/hazardlib/tests/calc/_conditioned_gmfs_test_data.py b/openquake/hazardlib/tests/calc/_conditioned_gmfs_test_data.py index 6bb7eb474bdd..e0f8aa5c0159 100644 --- a/openquake/hazardlib/tests/calc/_conditioned_gmfs_test_data.py +++ b/openquake/hazardlib/tests/calc/_conditioned_gmfs_test_data.py @@ -46,11 +46,8 @@ RUP = PointRupture( mag=MAG, - rake=0, tectonic_region_type=const.TRT.ACTIVE_SHALLOW_CRUST, hypocenter=Point(LON, LAT, DEP), - strike=0, - dip=90, occurrence_rate=1, temporal_occurrence_model=PoissonTOM(1.0), zbot=DEP, diff --git a/openquake/hazardlib/tests/gsim/mgmpe/modifiable_gmpe_test.py b/openquake/hazardlib/tests/gsim/mgmpe/modifiable_gmpe_test.py index cc353f8b2e98..db2b8a9624ba 100644 --- a/openquake/hazardlib/tests/gsim/mgmpe/modifiable_gmpe_test.py +++ b/openquake/hazardlib/tests/gsim/mgmpe/modifiable_gmpe_test.py @@ -31,11 +31,11 @@ class ModifiableGMPETest(unittest.TestCase): def test_AlAtik2015Sigma(self): + gmpe = valid.gsim('YenierAtkinson2015BSSA') params1 = {"tau_model": "global", "ergodic": False} params2 = {"tau_model": "cena", "ergodic": True} params3 = {} - gsims = [ModifiableGMPE(gmpe={'YenierAtkinson2015BSSA': {}}, - sigma_model_alatik2015=params) + gsims = [valid.modified_gsim(gmpe, sigma_model_alatik2015=params) for params in [params1, params2, params3]] cmaker = simple_cmaker(gsims, ['PGA']) ctx = cmaker.new_ctx(4) @@ -56,9 +56,10 @@ def test_AlAtik2015Sigma(self): aae(tau[2], 0.36855) # now test with_betw_ratio - cmaker.gsims[0] = ModifiableGMPE( + gsims[0] = ModifiableGMPE( gmpe={'Campbell2003': {}}, add_between_within_stds={'with_betw_ratio': 0.6}) + cmaker = simple_cmaker(gsims, ['PGA']) _mea, _sig, tau, phi = cmaker.get_mean_stds([ctx]) # (G,M,N) aae(tau[0], 0.44075136) aae(phi[0], 0.26445082) @@ -70,9 +71,9 @@ def test_AlAtik2015Sigma(self): def test_AkkarEtAlRjb2014(self): # check mean and stds - gsims = [ModifiableGMPE(gmpe={'AkkarEtAlRjb2014': {}}, - set_between_epsilon={'epsilon_tau': 0.5}), - valid.gsim('AkkarEtAlRjb2014')] + gmm = valid.gsim('AkkarEtAlRjb2014') + gsims = [valid.modified_gsim( + gmm, set_between_epsilon={'epsilon_tau': 0.5}), gmm] cmaker = simple_cmaker(gsims, ['PGA']) ctx = cmaker.new_ctx(4) ctx.mag = 6. @@ -102,9 +103,9 @@ def test_coefficients_as_dictionary(self): self.assertAlmostEqual(output_coeffs["XYZ"][SA(3.0)]["XYZ"], 3.0) def get_mean_stds(self, **kw): - gmpe_name = 'AkkarEtAlRjb2014' - gmm1 = ModifiableGMPE(gmpe={gmpe_name: {}}) - gmm2 = ModifiableGMPE(gmpe={gmpe_name: {}}, **kw) + gmm = valid.gsim('AkkarEtAlRjb2014') + gmm1 = valid.modified_gsim(gmm) + gmm2 = valid.modified_gsim(gmm, **kw) cmaker = simple_cmaker([gmm1, gmm2], ['PGA', 'SA(0.2)']) ctx = cmaker.new_ctx(4) ctx.mag = 6. @@ -183,9 +184,8 @@ def test(self): 'BaumontEtAl2018High2210IAVGDC30n7', 'FaccioliCauzzi2006']: - gmm = ModifiableGMPE(gmpe={gmpe_name: {}}, - apply_swiss_amplification={}) gmpe = valid.gsim(gmpe_name) + gmm = valid.modified_gsim(gmpe, apply_swiss_amplification={}) cmaker = simple_cmaker([gmm, gmpe], ['MMI']) ctx = self.get_ctx(cmaker) mea, _sig, _tau, phi = cmaker.get_mean_stds([ctx]) @@ -198,9 +198,8 @@ def test(self): 'EdwardsFah2013Foreland60Bars', 'ChiouYoungs2008SWISS01']: - gmm = ModifiableGMPE(gmpe={gmpe_name: {}}, - apply_swiss_amplification_sa={}) gmpe = valid.gsim(gmpe_name) + gmm = valid.modified_gsim(gmpe,apply_swiss_amplification_sa={}) cmaker = simple_cmaker([gmm, gmpe], ['SA(0.3)']) ctx = self.get_ctx(cmaker) mea, _sig, _tau, phi = cmaker.get_mean_stds([ctx]) diff --git a/openquake/qa_tests_data/scenario_risk/case_master/expected/aggrisk.csv b/openquake/qa_tests_data/scenario_risk/case_master/expected/aggrisk.csv index a25fbc88bf3f..6f89e8a88606 100644 --- a/openquake/qa_tests_data/scenario_risk/case_master/expected/aggrisk.csv +++ b/openquake/qa_tests_data/scenario_risk/case_master/expected/aggrisk.csv @@ -1,12 +1,12 @@ -#,,,"generated_by='OpenQuake engine 3.22.0-git4a826f5f28', start_date='2024-10-30T09:55:01', checksum=1711983385, investigation_time=None, risk_investigation_time=None" +#,,,"generated_by='OpenQuake engine 3.22.0-gite00dec3d72', start_date='2024-10-31T04:27:54', checksum=1711983385, investigation_time=None, risk_investigation_time=None" loss_type,rlz_id,loss_value,loss_ratio -business_interruption,0,1.80579E+03,1.28985E-01 -business_interruption,1,2.28556E+03,1.63254E-01 -contents,0,1.61214E+04,4.60610E-01 -contents,1,1.81703E+04,5.19151E-01 -nonstructural,0,3.01627E+04,2.87263E-01 -nonstructural,1,3.44602E+04,3.28192E-01 -structural,0,6.88028E+03,9.82897E-02 -structural,1,8.73611E+03,1.24802E-01 -occupants_night,0,5.37152E-02,1.27893E-03 -occupants_night,1,6.74443E-02,1.60582E-03 +business_interruption,0,1.79965E+03,1.28546E-01 +business_interruption,1,2.28696E+03,1.63355E-01 +contents,0,1.61509E+04,4.61454E-01 +contents,1,1.81652E+04,5.19007E-01 +nonstructural,0,3.01158E+04,2.86818E-01 +nonstructural,1,3.44231E+04,3.27839E-01 +structural,0,6.87515E+03,9.82165E-02 +structural,1,8.72907E+03,1.24701E-01 +occupants_night,0,5.37418E-02,1.27957E-03 +occupants_night,1,6.77008E-02,1.61192E-03 diff --git a/openquake/qa_tests_data/scenario_risk/case_master/expected/portfolio_loss.txt b/openquake/qa_tests_data/scenario_risk/case_master/expected/portfolio_loss.txt index eb5de88ed986..b0a70d975232 100644 --- a/openquake/qa_tests_data/scenario_risk/case_master/expected/portfolio_loss.txt +++ b/openquake/qa_tests_data/scenario_risk/case_master/expected/portfolio_loss.txt @@ -1,5 +1,5 @@ +------+-----------------------+----------+---------------+-----------+------------+ | loss | business_interruption | contents | nonstructural | occupants | structural | +------+-----------------------+----------+---------------+-----------+------------+ -| avg | 962.8 | 8_317 | 15_618 | 0.02857 | 3_672 | +| avg | 960.7 | 8_327 | 15_596 | 0.0286 | 3_669 | +------+-----------------------+----------+---------------+-----------+------------+ \ No newline at end of file diff --git a/openquake/qa_tests_data/scenario_risk/case_master/expected/portfolio_loss2.txt b/openquake/qa_tests_data/scenario_risk/case_master/expected/portfolio_loss2.txt index b1198e390881..368b5ad5ce16 100644 --- a/openquake/qa_tests_data/scenario_risk/case_master/expected/portfolio_loss2.txt +++ b/openquake/qa_tests_data/scenario_risk/case_master/expected/portfolio_loss2.txt @@ -1,5 +1,5 @@ +------+-----------------------+----------+---------------+-----------+------------+ | loss | business_interruption | contents | nonstructural | occupants | structural | +------+-----------------------+----------+---------------+-----------+------------+ -| avg | 7_035 | 34_291 | 63_202 | 0.2146 | 29_342 | +| avg | 7_015 | 34_304 | 62_873 | 0.2147 | 29_390 | +------+-----------------------+----------+---------------+-----------+------------+ \ No newline at end of file diff --git a/openquake/risklib/riskmodels.py b/openquake/risklib/riskmodels.py index dc1745596733..ac8606594f52 100644 --- a/openquake/risklib/riskmodels.py +++ b/openquake/risklib/riskmodels.py @@ -36,9 +36,8 @@ F32 = numpy.float32 F64 = numpy.float64 - -LTYPE_REGEX = '|'.join(lt for lt in scientific.LOSSTYPE - if '+' not in lt and '_ins' not in lt) +lts = numpy.concatenate([scientific.LOSSTYPE, scientific.PERILTYPE]) +LTYPE_REGEX = '|'.join(lt for lt in lts if '+' not in lt and '_ins' not in lt) RISK_TYPE_REGEX = re.compile(r'(%s)_([\w_]+)' % LTYPE_REGEX) @@ -56,7 +55,7 @@ def _assert_equal(d1, d2): def get_risk_files(inputs): """ :param inputs: a dictionary key -> path name - :returns: a pair (file_type, {risk_type: path}) + :returns: a dictionary "peril/kind/cost_type" -> fname """ rfs = {} job_ini = inputs['job_ini'] @@ -64,14 +63,22 @@ def get_risk_files(inputs): if key == 'fragility': # backward compatibily for .ini files with key fragility_file # instead of structural_fragility_file - rfs['fragility/structural'] = inputs[ + rfs['earthquake/fragility/structural'] = inputs[ 'structural_fragility'] = inputs[key] del inputs['fragility'] elif key.endswith(('_fragility', '_vulnerability', '_vulnerability_retrofitted')): match = RISK_TYPE_REGEX.match(key) if match: - rfs['%s/%s' % (match.group(2), match.group(1))] = inputs[key] - elif match is None: + kind = match.group(2) # fragility or vulnerability + value = inputs[key] + if isinstance(value, dict): # cost_type -> fname + peril = match.group(1) + for cost_type, fname in value.items(): + rfs[f'{peril}/{kind}/{cost_type}'] = fname + else: + cost_type = match.group(1) + rfs[f'earthquake/{kind}/{cost_type}'] = value + else: raise ValueError('Invalid key in %s: %s_file' % (job_ini, key)) return rfs @@ -96,29 +103,30 @@ def build_vf_node(vf): {'id': vf.id, 'dist': vf.distribution_name}, nodes=nodes) -def group_by_lt(funclist): +def group_by_peril(funclist): """ Converts a list of objects with attribute .loss_type into a dictionary - loss_type -> risk_function + peril -> loss_type -> risk_function """ - d = AccumDict(accum=[]) + ddic = AccumDict(accum=AccumDict(accum=[])) # peril -> lt -> rf for rf in funclist: - d[rf.loss_type].append(rf) - for lt, lst in d.items(): - if len(lst) == 1: - d[lt] = lst[0] - elif lst[1].kind == 'fragility': - # EventBasedDamageTestCase.test_case_11 - cf, ffl = lst - ffl.cf = cf - d[lt] = ffl - elif lst[1].kind == 'vulnerability_retrofitted': - vf, retro = lst - vf.retro = retro - d[lt] = vf - else: - raise RuntimeError(lst) - return d + ddic[rf.peril][rf.loss_type].append(rf) + for peril, dic in ddic.items(): + for lt, lst in dic.items(): + if len(lst) == 1: + dic[lt] = lst[0] + elif lst[1].kind == 'fragility': + # EventBasedDamageTestCase.test_case_11 + cf, ffl = lst + ffl.cf = cf + dic[lt] = ffl + elif lst[1].kind == 'vulnerability_retrofitted': + vf, retro = lst + vf.retro = retro + dic[lt] = vf + else: + raise RuntimeError(lst) + return ddic class RiskFuncList(list): @@ -132,7 +140,7 @@ def groupby_id(self): ddic = AccumDict(accum=[]) for rf in self: ddic[rf.id].append(rf) - return {riskid: group_by_lt(rfs) for riskid, rfs in ddic.items()} + return {riskid: group_by_peril(rfs) for riskid, rfs in ddic.items()} def get_risk_functions(oqparam): @@ -143,13 +151,13 @@ def get_risk_functions(oqparam): a list of risk functions """ job_ini = oqparam.inputs['job_ini'] - rmodels = AccumDict() + rmodels = AccumDict() # (peril, loss_type, kind) -> rmodel for key, fname in get_risk_files(oqparam.inputs).items(): - kind, loss_type = key.split('/') # ex. vulnerability/structural + peril, kind, loss_type = key.split('/') # ex. earthquake/vulnerability/structural rmodel = nrml.to_python(fname) if len(rmodel) == 0: raise InvalidFile(f'{job_ini}: {fname} is empty!') - rmodels[loss_type, kind] = rmodel + rmodels[peril, loss_type, kind] = rmodel if rmodel.lossCategory is None: # NRML 0.4 continue cost_type = str(rmodel.lossCategory) @@ -157,20 +165,16 @@ def get_risk_functions(oqparam): kind_ = kind.replace('_retrofitted', '') # strip retrofitted if not rmodel_kind.lower().startswith(kind_): raise ValueError( - 'Error in the file "%s_file=%s": is ' - 'of kind %s, expected %s' % ( - key, oqparam.inputs[key], rmodel_kind, - kind.capitalize() + 'Model')) + f'Error in the file "{key}_file={fname}": is ' + f'of kind {rmodel_kind}, expected {kind.capitalize() + "Model"}') if cost_type != loss_type: raise ValueError( - 'Error in the file "%s_file=%s": lossCategory is of ' - 'type "%s", expected "%s"' % - (key, oqparam.inputs[key], - rmodel.lossCategory, loss_type)) + f'Error in the file "{key}_file={fname}": lossCategory is of ' + f'type "{rmodel.lossCategory}", expected "{loss_type}"') cl_risk = oqparam.calculation_mode in ('classical', 'classical_risk') rlist = RiskFuncList() rlist.limit_states = [] - for (loss_type, kind), rm in sorted(rmodels.items()): + for (peril, loss_type, kind), rm in sorted(rmodels.items()): if kind == 'fragility': for (imt, riskid), ffl in sorted(rm.items()): if not rlist.limit_states: @@ -179,6 +183,7 @@ def get_risk_functions(oqparam): # limit states; this may change in the future assert rlist.limit_states == rm.limitStates, ( rlist.limit_states, rm.limitStates) + ffl.peril = peril ffl.loss_type = loss_type ffl.kind = kind rlist.append(ffl) @@ -187,6 +192,7 @@ def get_risk_functions(oqparam): # to make sure they are strictly increasing for (imt, riskid), rf in sorted(rm.items()): rf = rf.strictly_increasing() if cl_risk else rf + rf.peril = peril rf.loss_type = loss_type rf.kind = kind rlist.append(rf) @@ -212,12 +218,29 @@ def rescale(curves, values): return array +class PerilDict(dict): + """ + >>> pd = PerilDict({('earthquake', 'structural'): .23}) + >>> pd['structural'] + 0.23 + >>> pd['structurl'] + Traceback (most recent call last): + ... + KeyError: ('earthquake', 'structurl') + """ + def __getitem__(self, lt): + if isinstance(lt, tuple): + return dict.__getitem__(self, lt) + else: # assume lt is a loss_type string + return dict.__getitem__(self, ('earthquake', lt)) + + class RiskModel(object): """ Base class. Can be used in the tests as a mock. :param taxonomy: a taxonomy string - :param risk_functions: a dict (loss_type, kind) -> risk_function + :param risk_functions: a dict peril -> (loss_type, kind) -> risk_function """ time_event = None # used in scenario_risk compositemodel = None # set by get_crmodel @@ -232,11 +255,11 @@ def __init__(self, calcmode, taxonomy, risk_functions, **kw): if calcmode in ('classical', 'classical_risk'): self.loss_ratios = { lt: tuple(vf.mean_loss_ratios_with_steps(steps)) - for lt, vf in risk_functions.items()} + for lt, vf in risk_functions['earthquake'].items()} if calcmode == 'classical_bcr': self.loss_ratios_orig = {} self.loss_ratios_retro = {} - for lt, vf in risk_functions.items(): + for lt, vf in risk_functions['earthquake'].items(): self.loss_ratios_orig[lt] = tuple( vf.mean_loss_ratios_with_steps(steps)) self.loss_ratios_retro[lt] = tuple( @@ -244,7 +267,7 @@ def __init__(self, calcmode, taxonomy, risk_functions, **kw): # set imt_by_lt self.imt_by_lt = {} # dictionary loss_type -> imt - for lt, rf in self.risk_functions.items(): + for lt, rf in risk_functions['earthquake'].items(): if rf.kind in ('vulnerability', 'fragility'): self.imt_by_lt[lt] = rf.imt @@ -254,18 +277,21 @@ def loss_types(self): The list of loss types in the underlying vulnerability functions, in lexicographic order """ - return sorted(self.risk_functions) + return sorted(self.risk_functions['earthquake']) - def __call__(self, loss_type, assets, gmf_df, rndgen=None): + def __call__(self, assets, gmf_df, rndgen=None): meth = getattr(self, self.calcmode) - res = meth(loss_type, assets, gmf_df, rndgen) - return res # for event_based_risk this is a DataFrame (eid, aid, loss) + res = {(peril, lt): meth(peril, lt, assets, gmf_df, rndgen) + for peril in self.risk_functions for lt in self.loss_types} + # for event_based_risk this is a map loss_type -> DataFrame(eid, aid, loss) + return PerilDict(res) def __toh5__(self): return self.risk_functions, {'taxonomy': self.taxonomy} def __fromh5__(self, dic, attrs): vars(self).update(attrs) + assert 'earthquake' in dic, list(dic) self.risk_functions = dic def __repr__(self): @@ -273,7 +299,7 @@ def __repr__(self): # ######################## calculation methods ######################### # - def classical_risk(self, loss_type, assets, hazard_curve, rng=None): + def classical_risk(self, peril, loss_type, assets, hazard_curve, rng=None): """ :param str loss_type: the loss type considered @@ -288,7 +314,7 @@ def classical_risk(self, loss_type, assets, hazard_curve, rng=None): a composite array (loss, poe) of shape (A, C) """ n = len(assets) - vf = self.risk_functions[loss_type] + vf = self.risk_functions[peril][loss_type] lratios = self.loss_ratios[loss_type] imls = self.hazard_imtls[vf.imt] poes = hazard_curve[self.hazard_imtls(vf.imt)] @@ -302,7 +328,7 @@ def classical_risk(self, loss_type, assets, hazard_curve, rng=None): vf, imls, poes, lratios, self.investigation_time, rtime)] * n) return rescale(lrcurves, values) - def classical_bcr(self, loss_type, assets, hazard, rng=None): + def classical_bcr(self, peril, loss_type, assets, hazard, rng=None): """ :param loss_type: the loss type :param assets: a list of N assets of the same taxonomy @@ -315,7 +341,7 @@ def classical_bcr(self, loss_type, assets, hazard, rng=None): 'retrofitted is not defined for ' + loss_type) n = len(assets) self.assets = assets - vf = self.risk_functions[loss_type] + vf = self.risk_functions[peril][loss_type] imls = self.hazard_imtls[vf.imt] poes = hazard[self.hazard_imtls(vf.imt)] rtime = self.risk_investigation_time or self.investigation_time @@ -343,7 +369,7 @@ def classical_bcr(self, loss_type, assets, hazard, rng=None): for i, asset in enumerate(assets.to_records())] return list(zip(eal_original, eal_retrofitted, bcr_results)) - def classical_damage(self, loss_type, assets, hazard_curve, rng=None): + def classical_damage(self, peril, loss_type, assets, hazard_curve, rng=None): """ :param loss_type: the loss type :param assets: a list of N assets of the same taxonomy @@ -352,7 +378,7 @@ def classical_damage(self, loss_type, assets, hazard_curve, rng=None): where N is the number of points and D the number of damage states. """ - ffl = self.risk_functions[loss_type] + ffl = self.risk_functions[peril][loss_type] imls = self.hazard_imtls[ffl.imt] poes = hazard_curve[self.hazard_imtls(ffl.imt)] rtime = self.risk_investigation_time or self.investigation_time @@ -364,7 +390,7 @@ def classical_damage(self, loss_type, assets, hazard_curve, rng=None): for a in assets.to_records()]) return damages - def event_based_risk(self, loss_type, assets, gmf_df, rndgen): + def event_based_risk(self, peril, loss_type, assets, gmf_df, rndgen): """ :returns: a DataFrame with columns eid, eid, loss """ @@ -376,13 +402,13 @@ def event_based_risk(self, loss_type, assets, gmf_df, rndgen): else: val = assets['value-' + loss_type].to_numpy() asset_df = pandas.DataFrame(dict(aid=assets.index, val=val), sid) - vf = self.risk_functions[loss_type] + vf = self.risk_functions[peril][loss_type] return vf(asset_df, gmf_df, col, rndgen, self.minimum_asset_loss.get(loss_type, 0.)) scenario = ebrisk = scenario_risk = event_based_risk - def scenario_damage(self, loss_type, assets, gmf_df, rng=None): + def scenario_damage(self, peril, loss_type, assets, gmf_df, rng=None): """ :param loss_type: the loss type :param assets: a list of A assets of the same taxonomy @@ -396,7 +422,7 @@ def scenario_damage(self, loss_type, assets, gmf_df, rng=None): imt = self.imt_by_lt[loss_type] col = self.alias.get(imt, imt) gmvs = gmf_df[col].to_numpy() - ffs = self.risk_functions[loss_type] + ffs = self.risk_functions[peril][loss_type] damages = scientific.scenario_damage(ffs, gmvs).T return numpy.array([damages] * len(assets)) @@ -408,7 +434,7 @@ def scenario_damage(self, loss_type, assets, gmf_df, rng=None): # names of the parameter in the oqparam object. This is seen as a # feature, since it forces people to be consistent with the names, # in the spirit of the 'convention over configuration' philosophy -def get_riskmodel(taxonomy, oqparam, **extra): +def get_riskmodel(taxonomy, oqparam, risk_functions): """ Return an instance of the correct risk model class, depending on the attribute `calculation_mode` of the object `oqparam`. @@ -420,7 +446,7 @@ def get_riskmodel(taxonomy, oqparam, **extra): :param extra: extra parameters to pass to the RiskModel class """ - extra['hazard_imtls'] = oqparam.imtls + extra = {'hazard_imtls': oqparam.imtls} extra['investigation_time'] = oqparam.investigation_time extra['risk_investigation_time'] = oqparam.risk_investigation_time extra['lrem_steps_per_interval'] = oqparam.lrem_steps_per_interval @@ -430,7 +456,7 @@ def get_riskmodel(taxonomy, oqparam, **extra): if oqparam.calculation_mode == 'classical_bcr': extra['interest_rate'] = oqparam.interest_rate extra['asset_life_expectancy'] = oqparam.asset_life_expectancy - return RiskModel(oqparam.calculation_mode, taxonomy, **extra) + return RiskModel(oqparam.calculation_mode, taxonomy, risk_functions, **extra) # used only in riskmodels_test @@ -440,27 +466,28 @@ def get_riskcomputer(dic, alias): rc.asset_df = pandas.DataFrame(dic['asset_df']) rc.wdic = {} rfs = AccumDict(accum=[]) + steps = dic.get('lrem_steps_per_interval', 1) + mal = dic.get('minimum_asset_loss', {lt: 0. for lt in dic['loss_types']}) for rlk, func in dic['risk_functions'].items(): riskid, lt = rlk.split('#') rf = hdf5.json_to_obj(json.dumps(func)) if hasattr(rf, 'init'): rf.init() rf.loss_type = lt + rf.peril = 'earthquake' if getattr(rf, 'retro', False): rf.retro = hdf5.json_to_obj(json.dumps(rf.retro)) rf.retro.init() rf.retro.loss_type = lt rfs[riskid].append(rf) - steps = dic.get('lrem_steps_per_interval', 1) - mal = dic.get('minimum_asset_loss', {lt: 0. for lt in dic['loss_types']}) - for rlt, weight in dic['wdic'].items(): - riskid, lt = rlt.split('#') rm = RiskModel(dic['calculation_mode'], 'taxonomy', - group_by_lt(rfs[riskid]), + group_by_peril(rfs[riskid]), lrem_steps_per_interval=steps, minimum_asset_loss=mal) rm.alias = alias - rc[riskid, lt] = rm + rc[riskid] = rm + for rlt, weight in dic['wdic'].items(): + riskid, lt = rlt.split('#') rc.wdic[riskid, lt] = weight rc.loss_types = dic['loss_types'] rc.minimum_asset_loss = mal @@ -513,9 +540,13 @@ def read(cls, dstore, oqparam): if hasattr(dstore, 'get_attr'): # missing only in Aristotle mode, where dstore is an hdf5.File risklist.limit_states = dstore.get_attr('crm', 'limit_states') - df = dstore.read_df('crm', ['riskid', 'loss_type']) - for rf_json in df.riskfunc: + df = dstore.read_df('crm') + for i, rf_json in enumerate(df.riskfunc): rf = hdf5.json_to_obj(rf_json) + try: + rf.peril = df.loc[i].peril + except AttributeError: # in engine < 3.22 the peril was not stored + rf.peril = 'earthquake' lt = rf.loss_type if rf.kind == 'fragility': # rf is a FragilityFunctionList risklist.append(rf) @@ -644,8 +675,7 @@ def init(self): if oq.calculation_mode.endswith('_bcr'): # classical_bcr calculator for riskid, risk_functions in self.risklist.groupby_id().items(): - self._riskmodels[riskid] = get_riskmodel( - riskid, oq, risk_functions=risk_functions) + self._riskmodels[riskid] = get_riskmodel(riskid, oq, risk_functions) elif (any(rf.kind == 'fragility' for rf in self.risklist) or 'damage' in oq.calculation_mode): # classical_damage/scenario_damage calculator @@ -660,13 +690,11 @@ def init(self): self.damage_states = ['no_damage'] + list( self.risklist.limit_states) for riskid, ffs_by_lt in self.risklist.groupby_id().items(): - self._riskmodels[riskid] = get_riskmodel( - riskid, oq, risk_functions=ffs_by_lt) + self._riskmodels[riskid] = get_riskmodel(riskid, oq, ffs_by_lt) else: # classical, event based and scenario calculators for riskid, vfs in self.risklist.groupby_id().items(): - self._riskmodels[riskid] = get_riskmodel( - riskid, oq, risk_functions=vfs) + self._riskmodels[riskid] = get_riskmodel(riskid, oq, vfs) self.imtls = oq.imtls self.lti = {} # loss_type -> idx self.covs = 0 # number of coefficients of variation @@ -691,12 +719,14 @@ def init(self): if hasattr(rf, 'covs') and rf.covs.any(): self.covs += 1 self.curve_params = self.make_curve_params() + + # possibly set oq.minimum_intensity iml = collections.defaultdict(list) # ._riskmodels is empty if read from the hazard calculation for riskid, rm in self._riskmodels.items(): - for lt, rf in rm.risk_functions.items(): - if hasattr(rf, 'imt'): # vulnerability - iml[rf.imt].append(rf.imls[0]) + for lt, rf in rm.risk_functions['earthquake'].items(): + iml[rf.imt].append(rf.imls[0]) + if oq.aristotle: pass # don't set minimum_intensity elif sum(oq.minimum_intensity.values()) == 0 and iml: @@ -758,7 +788,7 @@ def make_curve_params(self): allratios = [] for taxo in sorted(self): rm = self[taxo] - rf = rm.risk_functions[loss_type] + rf = rm.risk_functions['earthquake'][loss_type] if loss_type in rm.loss_ratios: ratios = rm.loss_ratios[loss_type] allratios.append(ratios) @@ -850,12 +880,14 @@ def to_dframe(self): """ :returns: a DataFrame containing all risk functions """ - dic = {'riskid': [], 'loss_type': [], 'riskfunc': []} + dic = {'peril': [], 'riskid': [], 'loss_type': [], 'riskfunc': []} for riskid, rm in self._riskmodels.items(): - for lt, rf in rm.risk_functions.items(): - dic['riskid'].append(riskid) - dic['loss_type'].append(lt) - dic['riskfunc'].append(hdf5.obj_to_json(rf)) + for peril, rfdict in rm.risk_functions.items(): + for lt, rf in rfdict.items(): + dic['peril'].append(peril) + dic['riskid'].append(riskid) + dic['loss_type'].append(lt) + dic['riskfunc'].append(hdf5.obj_to_json(rf)) return pandas.DataFrame(dic) def __repr__(self): diff --git a/openquake/risklib/scientific.py b/openquake/risklib/scientific.py index a8a22a536c5a..fd162554ed17 100644 --- a/openquake/risklib/scientific.py +++ b/openquake/risklib/scientific.py @@ -45,6 +45,7 @@ 'losses', 'collapsed', 'injured', 'fatalities', 'homeless', 'non_operational'] +PERILTYPE = numpy.array(['earthquake', 'liquefaction', 'landslide']) LOSSTYPE = numpy.array('''\ business_interruption contents nonstructural structural occupants occupants_day occupants_night occupants_transit @@ -1663,12 +1664,12 @@ def __init__(self, crm, asset_df): self.wdic = {} tm = crm.tmap_df[crm.tmap_df.taxi == taxidx] country_str = getattr(asset_df, 'country', '?') - for lt in self.minimum_asset_loss: - for country, loss_type, riskid, weight in zip( - tm.country, tm.loss_type, tm.risk_id, tm.weight): + for country, loss_type, riskid, weight in zip( + tm.country, tm.loss_type, tm.risk_id, tm.weight): + for lt in self.minimum_asset_loss: if loss_type in ('*', lt): if country == '?' or country_str in country: - self[riskid, lt] = crm._riskmodels[riskid] + self[riskid] = crm._riskmodels[riskid] self.wdic[riskid, lt] = weight def output(self, haz, sec_losses=(), rndgen=None): @@ -1682,11 +1683,10 @@ def output(self, haz, sec_losses=(), rndgen=None): """ dic = collections.defaultdict(list) # lt -> outs weights = collections.defaultdict(list) # lt -> weights - for riskid, lt in self: - rm = self[riskid, lt] - out = rm(lt, self.asset_df, haz, rndgen) - weights[lt].append(self.wdic[riskid, lt]) - dic[lt].append(out) + for riskid, rm in self.items(): + for (peril, lt), res in rm(self.asset_df, haz, rndgen).items(): + weights[lt].append(self.wdic[riskid, lt]) + dic[lt].append(res) out = {} for lt in self.minimum_asset_loss: outs = dic[lt] @@ -1711,13 +1711,14 @@ def todict(self): """ rfdic = {} for rm in self.values(): - for lt, rf in rm.risk_functions.items(): - dic = ast.literal_eval(hdf5.obj_to_json(rf)) - if getattr(rf, 'retro', False): - retro = ast.literal_eval(hdf5.obj_to_json(rf.retro)) - dic['openquake.risklib.scientific.VulnerabilityFunction'][ - 'retro'] = retro - rfdic['%s#%s' % (rf.id, lt)] = dic + for peril, rfdict in rm.risk_functions.items(): + for lt, rf in rfdict.items(): + dic = ast.literal_eval(hdf5.obj_to_json(rf)) + if getattr(rf, 'retro', False): + retro = ast.literal_eval(hdf5.obj_to_json(rf.retro)) + dic['openquake.risklib.scientific.VulnerabilityFunction'][ + 'retro'] = retro + rfdic['%s#%s' % (rf.id, lt)] = dic df = self.asset_df dic = dict(asset_df={col: df[col].tolist() for col in df.columns}, risk_functions=rfdic, diff --git a/openquake/risklib/tests/riskmodels_test.py b/openquake/risklib/tests/riskmodels_test.py index 7e882a105666..e0991c7da7a7 100644 --- a/openquake/risklib/tests/riskmodels_test.py +++ b/openquake/risklib/tests/riskmodels_test.py @@ -205,6 +205,7 @@ def test1(self): 'RC#structural': {"openquake.risklib.scientific.VulnerabilityFunction": {"id": "RC", + "peril": 'earthquake', "loss_type": "structural", "imt": "PGA", "imls": [0.1, 0.2, 0.3, 0.5, 0.7], @@ -217,7 +218,8 @@ def test1(self): 'gmv_0': [.23, .31]} rc = riskmodels.get_riskcomputer(dic, alias={'PGA': 'gmv_0'}) print(toml.dumps(dic)) - self.assertEqual(dic, rc.todict()) + for k, v in rc.todict().items(): + self.assertEqual(dic[k], v) out = rc.output(pandas.DataFrame(gmfs)) print(out)