[cig-commits] [commit] master: Add mean slip/event rupture calc to QuakeLib (7630965)

[cig_noreply at geodynamics.org]

Repository : https://github.com/geodynamics/vq

On branch  : master
Link       : https://github.com/geodynamics/vq/compare/0ca7695073f8b223d815f10653613b02b6c47df8...197914f4e03e55d321cae45161deafb1a5ade706

>---------------------------------------------------------------

commit 76309653b8cec93f6fad92103d2a0f201457e28f
Author: Eric Heien <emheien at ucdavis.edu>
Date:   Wed Oct 29 16:18:15 2014 -0700

    Add mean slip/event rupture calc to QuakeLib


>---------------------------------------------------------------

76309653b8cec93f6fad92103d2a0f201457e28f
 examples/check_results.py | 30 ++++--------------------------
 quakelib/src/QuakeLibIO.h | 23 +++++++++++++++++++++++
 2 files changed, 27 insertions(+), 26 deletions(-)

diff --git a/examples/check_results.py b/examples/check_results.py
index 7e45a78..5235b3a 100755
--- a/examples/check_results.py
+++ b/examples/check_results.py
@@ -43,29 +43,11 @@ def check_self_consistent(events):
     return error
 
 def calc_mean_slip(events):
-    element_total_slips = {}
-    for event in events:
-        elements = event.getInvolvedElements()
-        element_sweep_slip_sums = {}
-        for elem in elements: element_sweep_slip_sums[elem] = 0
-        for sweep in event.getSweeps():
-            element_sweep_slip_sums[sweep._element_id] += sweep._slip
-            if not element_total_slips.has_key(sweep._element_id): element_total_slips[sweep._element_id] = 0.0
-            element_total_slips[sweep._element_id] += sweep._slip
-
-    mean_total_slip = 0.0
-    total_slips = [element_total_slips[elem_id] for elem_id in element_total_slips.keys()]
-
-    return sum(total_slips)/float(len(total_slips))
+    return sum([event.calcMeanSlip() for event in events])
 
 def calc_mean_interevent(events):
-    interevent_times = []
-    last_year = -1
-    for event in events:
-        if last_year > 0: interevent_times.append(event.getEventYear() - last_year)
-        last_year = event.getEventYear()
-
-    return sum(interevent_times)/float(len(interevent_times))
+    event_years = [event.getEventYear() for event in events]
+    return sum([event_years[i+1] - event_years[i] for i in range(len(event_years)-1)])/(len(event_years)-1)
 
 def calc_b_val(events):
     mags = [events.event_list[enum].magnitude for enum in events.event_list if events.event_list[enum].magnitude < 10 and events.event_list[enum].magnitude > 0]
@@ -83,11 +65,7 @@ def rupture_area_vs_mag(events):
     log_ra = []
     mag = []
     for event in events:
-        rupture_area = 0
-        element_areas = {}
-        for sweep in event.getSweeps():
-            element_areas[sweep._element_id] = sweep._area
-        rupture_area = sum([element_areas[bnum] for elem_num in element_areas])
+        rupture_area = event.calcEventRuptureArea()
 
         if not math.isnan(event.magnitude):
             log_ra.append(math.log10(rupture_area/1e6))
diff --git a/quakelib/src/QuakeLibIO.h b/quakelib/src/QuakeLibIO.h
index 1375af8..4397b70 100644
--- a/quakelib/src/QuakeLibIO.h
+++ b/quakelib/src/QuakeLibIO.h
@@ -915,6 +915,29 @@ namespace quakelib {
                 return _total_slip.size();
             };
 
+            //! Calculate the total slipped area in the event
+            double calcEventRuptureArea(void) const {
+                EventElementMap::const_iterator     it;
+                double                              rupture_area = 0;
+
+                for (it=_total_slip.begin(); it!=_total_slip.end(); ++it) rupture_area += it->second._area;
+
+                return rupture_area;
+            }
+
+            //! Calculate mean slip weighted by rupture area
+            double calcMeanSlip(void) const {
+                EventElementMap::const_iterator     it;
+                double                              sum_slip_area = 0, sum_area = 0;
+
+                for (it=_total_slip.begin(); it!=_total_slip.end(); ++it) {
+                    sum_slip_area += it->second._area*it->second._slip;
+                    sum_area += it->second._area;
+                }
+
+                return sum_slip_area/sum_area;
+            }
+
             void clear(void) {
                 _data._event_number = UNDEFINED_EVENT_ID;
                 _data._event_year = nan("");