jenndrei/bayhunter

McMC transdimensional Bayesian inversion of surface wave dispersion and receiver functions

bayhunter's Issues

numpy import fail

Hello, Jennifer
Thanks for this great tool. I was trying is on Ubuntu 64 bit machine in python2.7 environment and got the following error related to the numpy.
" ImportError: numpy.core.multiarray failed to import "

Could you please share the requirement file with versions of the python libs you are using.

Thanks
ajay

sudo python setup.py install reports an error

This is a great tool. I want to use it to do some valuable research, but now I have some small problems,I have installed numpy. Why keeps reminding me that there is no numpy

Thanks, Jennifer

Can you tell me how you solved this problem? Thank you very much

Thinning

Hi!
I just found this code when reading a paper, and I have one question: I read through the documentation, from what I understand every iteration is a perturbation of one of the parameters of the previous model. When doing this, you need to use some thinning (only put a model into the posterior every i iterations) else your successive models might be too strongly correlated. I didn't find any mention of thinning in the documentation, except a mention at the end of page 15:

every i-th chain model is saved

However, I can't find any other reference to this, does the code do thinning, is there a way to control the thinning parameter?
Regards,
Dorian Soergel

Suggestion for SynthObs.return_swddata()

First off, thank you for this code. It is great! We are using it for two projects right now and pulling out some very useful subsurface information.

We noticed that you don't save the any synthetic dispersion curves in the final c_*.npy files, so we are writing some tools to extract/compute the synthetic dispersion curves for a given model using SynthObs.return_swddata().

I have noticed two things that would be nice additions to this functions I think.

First, it would be good to check the inputs in return_swddata(). In particular I found that h (the depth axis) cannot be more than 100 elements. A useful error about this would be great.

Second, it would be nice to have another argument that allows the user to ask for a certain dispersion curves (e.g. Rayleigh wave group velocity). Right now, the code computes all 4 (phase and group for both Love and Rayleigh waves, which is a good default behavior).

I am happy to submit pull requests for these, but would like to know if they will actually get incorporated or not. Can you please let me know if this is something that you would incorporate?

Also, if anyone is still working on this code, it would be very useful to add more information about the inputs and outputs of the different methods. Right now, I feel that I am deciphering a lot of things on my own and this is slow. I know people will start using this code as more papers are published using it, but the user base would increase a lot if the code is a bit better documented. I am sure there is functionality that already exists in this code that I just don't know about, but because I don't know about it, I may be rewriting things from scratch. Regardless, thank you very much for putting this code online for us!!

A suggestion about freezing the model halfspace

Thanks to BayHunter, now one can apply the transdimensional MCMC inversion to receiver function and/or surface-wave dispersion conveniently. I find that in applying the inversion to a local-scale seismic data, it would be sometimes necessary to make the model halfspace consistent with the regional reference model. Thus I was thinking if it would be useful to add such an input, the halfspace parameters, to the inversion. Meanwhile, I come up with two ways to implement this:

1. In the get_vp_vs_h function in Models.py, one could modify the last component of the vp and vs arrays. By doing this, one modifies the model halfspace directly during the inversion.
2. One could also add the halfspace parameters to the thkm, vpm, vsm, and rhom arrays in the get_modelvectors function of surf96_modsw.py. One also needs to add 1 to nlayer in run_model. This only works for surface-wave dispersion and actually adds a layer to the inversion model. Thus I do not recommend this way.

I am not sure if these implementing would be added. I wrote these down in case someone will need to freeze the model halfspace in BayHunter.

Missing import os in SynthObs.py

I was trying to save synthetic dispersion curves to a given file name. I got the following error.

SynthObs.save_data(syn_disp, outfile=disp_file)

File "/Users/dmikesell/anaconda3/envs/tomo/lib/python3.7/site-packages/BayHunter/SynthObs.py", line 108, in save_data
name, ext = os.path.splitext(outfile)
NameError: name 'os' is not defined

I will add os to the import list in this file. I will also fix an error in the "ext" and the ".". When providing a file name, it ends up having two dots in the extension because "ext" already contains a dot.

Not able to use the package for 50 mtr depth.

Dear Jennifer

I am trying to use Bayhunter package to invert some near seismic data for SWD only. The maximum depth of the data that i have is 50 m. There are some issues in the plotting package with respect to binning size (default value for which is 2kms. in the package)

There may be some other issues in relation to my work, i seek your help regarding the same.
Thank you very much.

install rfmini

If anybody else stumbles upon some errors during installation of rfmini with recent python versions, here is how I solved mine.

I needed to delete the file rfmini.c. Furthermore I edited the setup.py script of rfmini (line with include_dirs):

from distutils.core import setup
from distutils.extension import Extension
from Cython.Distutils import build_ext
import numpy as np

ext_modules=[
    Extension("rfmini", ["rfmini.pyx","greens.cpp","model.cpp","pd.cpp","synrf.cpp","wrap.cpp","fork.cpp"],
                  include_dirs=[np.get_include()])
]

setup(
  name = 'RfMini',
  cmdclass = {'build_ext': build_ext},
  ext_modules = ext_modules,
)

In my conda environment I installed the following packages:

conda install cython cxx-compiler

After that installation with pip was successful.

Bug in synthetic dispersion curve save_data() method

We have been using BayHunter for a near-surface study (<100m). We have been using the built-in routine to SynthObs.return_swddata() to compute synthetic dispersion curves for the mean and median models. With one of our data so far, we have run into a tricky little bug. For the mean Vs model, we can compute synthetic curve with no problem; however, we get an error with the median model. This occurs for the Love wave. The model is a thin fast layer over a halfspace. The error is strange because it does not exist on a mac, but does occur on our linux cluster.

$ python write_npy2csv.py 
mean
dict_keys(['rdispph', 'rdispgr', 'ldispph', 'ldispgr'])
ref = rdispph, type(x) = <class 'numpy.ndarray'>
ref = rdispgr, type(x) = <class 'numpy.ndarray'>
ref = ldispph, type(x) = <class 'numpy.ndarray'>
ref = ldispgr, type(x) = <class 'numpy.ndarray'>
median
dict_keys(['rdispph', 'rdispgr', 'ldispph', 'ldispgr'])
ref = rdispph, type(x) = <class 'numpy.ndarray'>
ref = rdispgr, type(x) = <class 'numpy.ndarray'>
ref = ldispph, type(x) = <class 'numpy.float64'>
Traceback (most recent call last):
  File "write_npy2csv.py", line 120, in <module>
    extract_data(depth_int, syn_per, basepath, plot=True)
  File "write_npy2csv.py", line 57, in extract_data
    SynthObs.save_data(syn_disp, outfile=disp_file)
  File "/home/dylanmikesell/anaconda3/lib/python3.7/site-packages/BayHunter/SynthObs.py", line 116, in save_data
    if len(x) > 0: # check if targets exist
TypeError: object of type 'numpy.float64' has no len()

You can see in the last print statement that the ldispph in save_data() is all of a sudden a float64 rather than an ndarray. This does not happen on the mac...type(x) remains an ndarray for this model.

Are there any ideas of what could cause this behavior?