MLh calculation for intermediate-depth earthquakes

[kpapaza]

Dear friends from SED,

we (and of course alomost everybody in Seismology in Europe) are using Seiscomp, an as part of our standard process we would like to be able to calculate local magnitude (MLh) even for depths greater than 80 km. This issue is critical for countries like Greece, Italy, Spain, Romania, etc., that have intermediate-depth seismicity, with events with depths > 80km

Unfortunately, this limit (80km) is hard coded in a couple of files (Mlh.cpp, ml.cpp), which means that we have to modify them and recompile Seiscomp in order to solve the problem. With SC3 we used to follow this procedure, but now that we’ re migrating to SC5 we’d like to avoid it as much as possible. Re-compiling creates several problems, especially regarding patches (as well as other issues), that we would like to avoid. I have contacted several othe people (from NOA-Athens, NIEP-Romania) and they only solved the problem by a) Changing the code and, b) Recompiling Seiscomp. Jans Becker from Seiscomp has confirmed that this is the only option (currently)

Since MLh is part of the SED contributions and a plugin, we would like to ask you if it woul;d be possible to either relax this constraint (remove completelya depth limitation), or allow the user to modify a relevant parameter in a .cfg file, either for the existing Seiscomp5 (in a future update), or in the new Seiscomp6.

Thank you in advance for your help and understanding

Costas Papazachos from Aristotle Univ. Thessaloniki network (HT)

PS. I have not commented on the science issue of the validity and use of MLh for intermediate-depth events, this is a different story...

[luca-s]

Hi @kpapaza . I don't mind make the parameter configurable, but I need to double check with my colleagues, the experts of MLh, to verify if it makes sense at all to compute MLh at those depths. MLh was designed for Switzerland and so I would like to make sure it is scientifically sensible to allow such additional configuration.

[kpapaza]

Dear Luca,

    thank you for your email and your willingness to help with the MLh issue. Being in a country with intermediate-depth events, I can already tell you that the significance of MLh for intermediate-depth events is really very poorly justified (from a scientific point of view).

I am in the lucky position to have worked as a researcher with these types of events with my students for many years, and I can tell you several practical reasons:

a) These events have tremendous along-arc/back-arc variability (one order of magnitude smaller amplitudes at back-arc stations for the Aegean arc deep ebents, for events with h>100km). While this was known since the first time my father addressed this issue (Papazachos and Comninakis, 1971, https://agupubs.onlinelibrary.wiley.com/doi/pdf/10.1029/JB076i035p08517?casa_token=i3W5KUQtLmoAAAAA:CReFkEJE28lThQ_e5wjO23eBLUF7JMQcjQdVvwfiPzUE4yjxuffjsL6rjRzrxGy920COiCfR0i62CpCL), we demonstrated this also quantitatively with my student Andreas Skarlatoudis using strong and weak motion data (e.g. Skarlatoudis et al., 2013, https://pubs.geoscienceworld.org/ssa/bssa/article-abstract/103/3/1952/331704/Ground-Motion-Prediction-Equations-of-Intermediate). This strong variability due to the mantle wedge attenuation, which has a very low Qp/Qs (e.g. see Ventouzi et al., https://academic.oup.com/gji/article/215/1/635/5059579), also has strong implications for their spectral content, as seen from stochastic modelling (e.g. see Kkallas et al., 2018, https://pubs.geoscienceworld.org/ssa/bssa/article/108/2/946/525973/Stochastic-Strong-Ground-Motion-Simulation-of-the?casa_token=vjAlp8begXYAAAAA%3aEt80q3a7skpPgz8JarUcuOUZEPuCq3KhjRtssXdYo86hZ4KaeRGheweCyxtlRhzsPRf_gWi7). Finally, this is also clearly seen in historical and recent intermediate-depth damage data modeling (e.g., Kkallas et al., 2018, https://link.springer.com/article/10.1007/s10518-018-0342-8). Sorry for bombing you with the publications from my group, I simply want to let you know that we have worked extensively with such events and are very aware of the problems and challenges that their processing poses.

The main problem is that the S-wavetrain records fro such events are "killed" in the back-arc area, while being amplified (channeled along the slab) for the along-arc stations, resulting in MLh differences of the order of 1 magnitude unit. Not to mention that standard MLh corrections (logA0), are built for shallow events, and do not make much sense for such deep events.

b) If this is the situation, why do we insist on computing MLh? There are 2 main reasons:

1) Computing robust magnitudes for deep events is a challenge (see our work with A. Tsambas for global events using mb, Ms and MJMA against Mw, https://pubs.geoscienceworld.org/ssa/bssa/article-abstract/106/2/418/332150/Global-Magnitude-Scaling-Relations-for).

2) In Greece we have computed since the 1950s a large number of ML from horizontal components for deep events. Being able to keep computing MLh today gives us a "time-continuation", which will allow as (at some point) to apply proper statistical corrections, and re-compute equivalent moment magnitudes for older events, which had been recorded on paper and for which nothing can be done nowdays, in other words to generate a homogeneous catalogue.

I hope this helps your colleagues to decide about this issue.

Best regards

Costas Papazachos

On 2023-03-15 18:06, luca-s wrote:

Hi @kpapaza https://github.com/kpapaza . I don't mind make the parameter configurable, but I need to double check with my colleagues, the experts of MLh, to verify if it makes sense at all to compute MLh at those depths. MLh was designed for Switzerland and so I would like to make sure it is scientifically sensible to allow such additional configuration.

— Reply to this email directly, view it on GitHub https://github.com/swiss-seismological-service/sed-SeisComP-contributions/issues/15#issuecomment-1470318293, or unsubscribe https://github.com/notifications/unsubscribe-auth/A6QGJLUCRE2MLJKIB55MG6DW4HSHFANCNFSM6AAAAAAV3XFXIM. You are receiving this because you were mentioned.Message ID: @.***>

-- Costas Papazachos Professor of Geophysics, Geophysical Laboratory, Aristotle University of Thessaloniki, PO Box 352-1, GR-54124, Thessaloniki, GREECE @.*** Tel.: + 30 2310998510, Fax: + 30 2310998528

(Greek Web Page)http://geophysics.geo.auth.gr/new_web_site_2007/GR/dep/DEP_PERSONAL/costas.html (English Web Page)http://geophysics.geo.auth.gr/new_web_site_2007/EN/dep/DEP_PERSONAL/costas.html

Email Disclaimer The information in this email is confidential and is intended solely for the addressee(s). If you have received this transmission in error, and you are not an intended recipient, be aware that any disclosure, copying, distribution or use of this transmission or its contents is prohibited. Furthermore, you are kindly requested to send us back the original message to the sender’s address and delete the message from your system immediately. Internet communications are not secure and therefore AUTH does not accept legal responsibility for the contents of this message and for any damage whatsoever that is caused by viruses being passed. Thank You, Aristotle University of Thessaloniki

Διευκρίνιση ηλεκτρονικού ταχυδρομείου Οι πληροφορίες που συμπεριλαμβάνονται σε αυτό το μήνυμα είναι εμπιστευτικές και η χρήση τους επιτρέπεται μόνον από τον αναφερόμενο παραλήπτη. Εάν έχετε λάβει το παρόν μήνυμα από λάθος και δεν είστε ο προοριζόμενος παραλήπτης, σας ενημερώνουμε ότι αποκάλυψη, αναπαραγωγή, διανομή ή οποιασδήποτε άλλης μορφής χρήση των περιεχομένων του παρόντος μηνύματος απαγορεύεται. Επίσης παρακαλείσθε να αποστείλετε το αρχικό μήνυμα στην διεύθυνση του αποστολέα, καθώς και στη συνέχεια να διαγράψετε το μήνυμα από το σύστημά σας. Η επικοινωνία μέσω Internet δεν είναι ασφαλής και επομένως το ΑΠΘ δεν φέρει ευθύνη για οποιαδήποτε θετική ή αποθετική ζημιά που προκλήθηκε από την χρήση του παρόντος ή των συνημμένων του λόγω ιών που έχουν περάσει σε αυτά. Σας Ευχαριστούμε, Αριστοτέλειο Πανεπιστήμιο Θεσσαλονίκης

[luca-s]

Thank you for the detailed explanation @kpapaza . Knowing that you have a good understanding of the implications makes me more comfortable in allowing the configuration of this parameter. It would probably makes just sense to add a warning in the new parameter description stating something on the line "Change the default value of 80km only if you know what you are doing" ;)

Since you have already applied the changes on the code, would you like to provide a pull request? Or did you simply modified the hardcoded values?

[kpapaza]

Dear Luca,

    thank you for your kind email. A warning will really be helpful!

We simply have modified the existing hardcoded values. I will ask our sysasmin, Petros Triantafyllidis, and Seiscomp admin, Odysseas Galanis, if feel comfortable to modify the code properly and submit a pull request and let you know.

Best regards

Costas

On 2023-03-16 10:58, luca-s wrote:

Thank you for the detailed explanation @kpapaza https://github.com/kpapaza . Knowing that you have a good understanding of the implications makes me more comfortable in allowing the configuration of this parameter. It would probably makes just sense to add a warning in the new parameter description stating something on the line "Change the default value of 80km only if you know what you are doing" ;)

Since you have already applied the changes on the code, would you like to provide a pull request? Or did you simply modified the hardcoded values?

— Reply to this email directly, view it on GitHub https://github.com/swiss-seismological-service/sed-SeisComP-contributions/issues/15#issuecomment-1471551451, or unsubscribe https://github.com/notifications/unsubscribe-auth/A6QGJLQ3JFKRW36Y5HA25C3W4LI4ZANCNFSM6AAAAAAV3XFXIM. You are receiving this because you were mentioned.Message ID: @.***>

-- Costas Papazachos Professor of Geophysics, Geophysical Laboratory, Aristotle University of Thessaloniki, PO Box 352-1, GR-54124, Thessaloniki, GREECE @.*** Tel.: + 30 2310998510, Fax: + 30 2310998528

(Greek Web Page)http://geophysics.geo.auth.gr/new_web_site_2007/GR/dep/DEP_PERSONAL/costas.html (English Web Page)http://geophysics.geo.auth.gr/new_web_site_2007/EN/dep/DEP_PERSONAL/costas.html

Email Disclaimer The information in this email is confidential and is intended solely for the addressee(s). If you have received this transmission in error, and you are not an intended recipient, be aware that any disclosure, copying, distribution or use of this transmission or its contents is prohibited. Furthermore, you are kindly requested to send us back the original message to the sender’s address and delete the message from your system immediately. Internet communications are not secure and therefore AUTH does not accept legal responsibility for the contents of this message and for any damage whatsoever that is caused by viruses being passed. Thank You, Aristotle University of Thessaloniki

Διευκρίνιση ηλεκτρονικού ταχυδρομείου Οι πληροφορίες που συμπεριλαμβάνονται σε αυτό το μήνυμα είναι εμπιστευτικές και η χρήση τους επιτρέπεται μόνον από τον αναφερόμενο παραλήπτη. Εάν έχετε λάβει το παρόν μήνυμα από λάθος και δεν είστε ο προοριζόμενος παραλήπτης, σας ενημερώνουμε ότι αποκάλυψη, αναπαραγωγή, διανομή ή οποιασδήποτε άλλης μορφής χρήση των περιεχομένων του παρόντος μηνύματος απαγορεύεται. Επίσης παρακαλείσθε να αποστείλετε το αρχικό μήνυμα στην διεύθυνση του αποστολέα, καθώς και στη συνέχεια να διαγράψετε το μήνυμα από το σύστημά σας. Η επικοινωνία μέσω Internet δεν είναι ασφαλής και επομένως το ΑΠΘ δεν φέρει ευθύνη για οποιαδήποτε θετική ή αποθετική ζημιά που προκλήθηκε από την χρήση του παρόντος ή των συνημμένων του λόγω ιών που έχουν περάσει σε αυτά. Σας Ευχαριστούμε, Αριστοτέλειο Πανεπιστήμιο Θεσσαλονίκης

[luca-s]

Dear @kpapaza ,

to avoid that the next seiscomp release ships without this feature I moved on and implemented it myself. I thought it was easier for me to implement it since I always work on Seiscomp code.

I would like to ask you to double check that my changes work fine in your system and that you are getting the same MLh values as you currently do. You can find my change is in MLhmaxDepth branch. if you confirm it works fine I will merge it to master and it will be come available in the next SeisComP release.

The MLh magnitude bindings have now a maxDepth parameter.

However that is not enough, also the amplitude MLh bindings have to be configured with the desired maxDepth. This can be achieved by adding an "Amplitude type profile" with the name MLh (the name must be the same as the amplitude type). In there you'll find the maxDepth parameter for the amplitude.

Please let me know if everything works as expected .

Best Luca

[gempa-jabe]

Generally you could also easily implement that in a locale (localization of magnitudes). Each locale receives a minimum depth and maximum depth and a region (and some more). If a locale is set, then use its values to compute the magnitude. Your changes will also work, no problem, locales might be a more elegant solution although pretty new to SeisComP and not widely known yet.

[luca-s]

@gempa-jabe your proposal seems much more elegant. As you said, "It is not widely known yet", and indeed I was not aware of it :)

[gempa-jabe]

The cool thing about locales is that you can also bind additional parameters for your particular magnitude and then customize the attenuation table as well. More fancy stuff with magnitudes to come ...

[kpapaza]

Dear @luca-s, thank you for the code modification and the introduction of this new feature. We will test it and let you know if we face any issues or changes, in comparison to the previous version and its calculations.

Dear @gempa-jabe, thank you for your input. I was really not aware that it is possible to use such a feature (locale) to introduce a custom magnitude with its own attenuation table. This is not only elegant but also extremely useful, if e.g. localized attenuation tables have been proposed fot an area, country, region, etc.

Costas

[luca-s]

@kpapaza just to be clear, the changes are on a separate branch and until I merge them on the master branch of this code repository they will not appear in the next release. For this reason I would like to ask you to confirm that the changes in that branch work as expected for you. If you verify that I would then merge them to master to make them available in the next SeisComP release. Would you be able to do this test?

[kpapaza]

@luca-s Thank you for the clarification, we will test this at our test Seiscomp installation, check how it performs using historical data (older deep events, with h > 80km) and get back to you when we have something solid...

[ogalanis]

@luca-s I installed your branch and calculated some magnitudes, preserving the same configuration for everything else. The results are, as far as I can tell, and as expected, identical to our modified installation. I would recommend you merge it to the master branch. Thanks again, on behalf of our network and other users with the same issue.

cc: @kpapaza

[luca-s]

Thank you for reporting back @ogalanis . I'll move on and merge this feature then. Unfortunately I see that SeicComP 5.4.0 has been already tagged so you'll see this feature on the next release.

[ogalanis]

@luca-s Hi. I hate to resurrect this issue, but it seems that somewhere along the way (certainly after 5.5.6) the feature broke. Specifically, I observed that in release 6.4.3 (probably in earlier releases as well) when magnitudes.MLh.maxDepth and amplitudes.MLh.maxDepth are set to >= 80, it crashes scmag. Amplitudes seem to be produced correctly without any problem, but as soon as scmag gets a new origin and tries to calculate MLh station magnitudes it just crashes silently. Is it an issue with SED's contribution, the main SC (more likely) or some incompatibility between the two? If I forgot to report anything that might help, like log files, please ask me.

[luca-s]

Hi @ogalanis, let's find out where the code crashes. The easiest way is gdb.

Stop scmag

seiscomp stop scmag

Start scmag via gdb

gdb --args scmag

In the gdb session start the module (r):

Reading symbols from scmag...
(gdb) r
Starting program: scmag

When the module crashes, type the backtrace gdb command (bt) to see where it crashed and report it here:

(gdb) bt

[ogalanis]

Thanks @luca-s . I really appreciate your help and your time. Here's what happened: I set the maxDepth to 90 for just some stations, and sure enough it crashed. Here is the backtrace. I tried a few times to see if it is consistent:

#0  __pthread_kill_implementation (threadid=<optimized out>, signo=signo@entry=6,
no_tid=no_tid@entry=0) at pthread_kill.c:44
#1  0x00007ffff4e8b9b3 in __pthread_kill_internal (signo=6, threadid=<optimized out>)
at pthread_kill.c:78
#2  0x00007ffff4e3e646 in __GI_raise (sig=sig@entry=6) at ../sysdeps/posix/raise.c:26
#3  0x00007ffff4e287f3 in __GI_abort () at abort.c:79
#4  0x00007ffff52a1b21 in __gnu_cxx::__verbose_terminate_handler ()
at ../../../../libstdc++-v3/libsupc++/vterminate.cc:95
#5  0x00007ffff52ad52c in __cxxabiv1::__terminate (handler=<optimized out>)
at ../../../../libstdc++-v3/libsupc++/eh_terminate.cc:48
#6  0x00007ffff52ad597 in std::terminate ()
at ../../../../libstdc++-v3/libsupc++/eh_terminate.cc:58
#7  0x00007ffff52ad7f9 in __cxxabiv1::__cxa_throw (obj=<optimized out>,
tinfo=0x7ffff57fabf0 <typeinfo for fmt::v9::format_error>,
dest=0x7ffff57ba3fe fmt::v9::format_error::~format_error())
at ../../../../libstdc++-v3/libsupc++/eh_throw.cc:95
#8  0x00007ffff57ba042 in fmt::v9::detail::throw_format_error(char const*) ()
from /home/sysop/seiscomp/lib/libseiscomp_config.so
#9  0x00007ffff6a7f188 in fmt::v9::appender fmt::v9::detail::write_int_noinline<char, fmt::v9::appender, unsigned int>(fmt::v9::appender, fmt::v9::detail::write_int_arg<unsigned int>, fmt::v9::basic_format_specs<char> const&, fmt::v9::detail::locale_ref) ()
from /home/sysop/seiscomp/lib/libseiscomp_core.so.16
#10 0x00007ffff6a78048 in fmt::v9::appender fmt::v9::detail::printf_arg_formatter<fmt::v9::appender, char>::operator()<unsigned int, 0>(unsigned int) ()
from /home/sysop/seiscomp/lib/libseiscomp_core.so.16
#11 0x00007ffff6a74f56 in void fmt::v9::detail::vprintf<char, fmt::v9::basic_printf_context<fmt::v9::appender, char> >(fmt::v9::detail::buffer<char>&, fmt::v9::basic_string_view<char>, fmt::v9::basic_format_args<fmt::v9::basic_printf_context<fmt::v9::appender, char> >) () from /home/sysop/seiscomp/lib/libseiscomp_core.so.16
#12 0x00007ffff6a72f0f in std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > fmt::v9::vsprintf<fmt::v9::basic_string_view<char>, char>(fmt::v9::basic_string_view<char> const&, fmt::v9::basic_format_args<fmt::v9::basic_printf_context<std::conditional<std::is_same<fmt::v9::type_identity<char>::type, char>::value, fmt::v9::appender, std::back_insert_iterator<fmt::v9::detail::buffer<fmt::v9::type_identity<char>::type> > >::type, fmt::v9::type_identity<char>::type> >) ()
from /home/sysop/seiscomp/lib/libseiscomp_core.so.16
#13 0x00007ffff6a6f64f in Seiscomp::Logging::VPublish(Seiscomp::Logging::PublishLoc*, Seiscomp::Logging::Channel*, fmt::v9::basic_string_view<char>, fmt::v9::basic_format_args<fmt::v9::basic_printf_context<fmt::v9::appender, char> >) ()
from /home/sysop/seiscomp/lib/libseiscomp_core.so.16
--Type <RET> for more, q to quit, c to continue without paging--
#14 0x00007ffff511defa in void Seiscomp::Logging::PublishP<char [132], double&, int>(Seiscomp::Logging::PublishLoc*, Seiscomp::Logging::Channel*, char const (&) [132], double&, int&&) () from /home/sysop/seiscomp/share/plugins/mlh.so
#15 0x00007ffff5119a2a in (anonymous namespace)::MagnitudeProcessor_ML::setup(Seiscomp::Processing::Settings const&) () from /home/sysop/seiscomp/share/plugins/mlh.so
#16 0x000000000044431c in Seiscomp::Magnitudes::MagTool::computeStationMagnitude(Seiscomp::DataModel::Amplitude const*, Seiscomp::DataModel::Origin const*, Seiscomp::DataModel::SensorLocation const*, double, double, std::vector<Seiscomp::Magnitudes::MagTool::MagnitudeEntry, std::allocatorSeiscomp::Magnitudes::MagTool::MagnitudeEntry >&) ()
#17 0x000000000044dd46 in Seiscomp::Magnitudes::MagTool::feed(Seiscomp::DataModel::Amplitude*, bool, bool) ()
#18 0x000000000046f316 in MagToolApp::addObject(std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, Seiscomp::DataModel::Object*) ()
#19 0x00007ffff7c9e062 in Seiscomp::Client::Application::handleNotifier(Seiscomp::DataModel::Notifier*) () from /home/sysop/seiscomp/lib/libseiscomp_client.so.16
#20 0x00007ffff7c9dede in Seiscomp::Client::Application::handleMessage(Seiscomp::Core::Message*) () from /home/sysop/seiscomp/lib/libseiscomp_client.so.16
#21 0x000000000046f0f9 in MagToolApp::handleMessage(Seiscomp::Core::Message*) ()
#22 0x00007ffff7c962a8 in Seiscomp::Client::Application::dispatch(Seiscomp::Core::BaseObject*) () from /home/sysop/seiscomp/lib/libseiscomp_client.so.16
#23 0x00007ffff7c95515 in Seiscomp::Client::Application::processEvent() ()
from /home/sysop/seiscomp/lib/libseiscomp_client.so.16
#24 0x00007ffff7c95300 in Seiscomp::Client::Application::run() ()
from /home/sysop/seiscomp/lib/libseiscomp_client.so.16
#25 0x000000000046eef8 in MagToolApp::run() ()
#26 0x00007ffff6af6b25 in Seiscomp::System::Application::exec() ()
from /home/sysop/seiscomp/lib/libseiscomp_core.so.16
#27 0x000000000046a730 in main ()

[luca-s]

It could be this line. The logging API changed in SeisComP 6.

Could you try recompiling the code after replacing DEPTH_MAX with double(DEPTH_MAX) at that line?

[ogalanis]

This did the trick, indeed! However, we are using the pre-compiled binaries for CentOS for our production systems. Can this line be changed in the next release?

[luca-s]

Sure, we certainly have to push the fix to master and made it available to the next release.

@gempa-jabe I committed the fix on master. Are you planning a bug fix release any time soon?

[gempa-jabe]

I have checked the logging issue and it is throwing an exception with "invalid format specifier" because it expects "%d". To solve the issue it would help to either explicitly define MAX_DEPTH as double or correct the output specifier.

#define DEPTH_MAX 80.0

Your fix works as well but is unnecessarily complex (casting and double formatting) and maybe more difficult to understand later on.

We are preparing for the 6.5 release which takes some more time due to the holiday season. Of course your fix will be part of it.

[luca-s]

@gempa-jabe I understand the cause of the crash and I initially thought of fixing it like your suggestion, but I didn't want to leave the code open to other crashes in case somebody changed the value of DEPTH_MAX again. However feel free to change the code as you prefer, I don't have a strong opinion on that. Thank you.

[gempa-jabe]

Thanks for the reminder of the "external" change. Given that use-case, it makes sense to keep it like it is.