geophysics/mtpy

collection of python tools for standard MT data processing

MTpy

A Python toolbox for using magnetotelluric data.

Note: This project has moved to https://github.com/MTgeophysics/mtpy

Installation

$ python setup.py install

If you have setuptools installed, we suggest you use:

$ python setup.py develop

Documentation

See http://geophysics.github.com/mtpy

License

MTpy is licensed under the GPL version 3

The license agreement is contained in the repository and should be kept together with the code.

Copyright: Lars Krieger, Jared Peacock, Kent Inverarity, Alison Kirkby; Adelaide, 2013

Conventions

1. MTpy uses E- and B-fields (although the sensors may be confusingly named as H-sensors in EDI files)
2. [E] = microvolts/meter (muV/m)
3. [B] = nanotesla (nT)
4. [Z] = [E]/[B] = km/s
5. Apparent resistivty rho = 0.2 * T * |Z|^2 (in Ohm m)
6. Angles are given in degrees (mod 360)
7. EDI files can contain data in Z- or rho/phi-form
8. EDI files contain data from one station only
9. Coordinates are handled in decimal degrees (converted when reading)
10. Time stamps refer to UTC
11. Internal coordinates: X = North, Y = East
12. Rotations are interpreted clockwise (mathematically negative)
13. 0 degrees azimuth = North

mtpy's Issues

setuptools' 'entry_points' cannot call a module as a script

> occam2d_gui.exe
Traceback (most recent call last):
  File "C:\Python27\Scripts\occam2d_gui-script.py", line 9, in <module>
    load_entry_point('MTpy==0.0.0', 'console_scripts', 'occam2d_gui')()
TypeError: 'module' object is not callable

The fix is to use the convention for script/modules:

import sys

def main(args):
    # do stuff

if __name__ == '__main__':
    main(sys.argv[1:])

http://effbot.org/pyfaq/tutor-what-is-if-name-main-for.htm
http://stackoverflow.com/a/419189

I'll change the existing scripts to do this.

Ploting phase tensor ellipses over a geological map

Dear All

I have a geological map of my study area and I want to plot Phase tensor ellipse map
over the background of geological map. Please suggest me necessary steps.

Best Regards

there is a bug in ws3dinv.py

import mtpy.modeling.ws3dinv as ws
import os
.............................
mm = ws.WS3DModelManipulator(initial_fn=wsmesh.initial_fn)

here is no WS3DModelManipulator function in the module of ws3dinv.py and it will get a error. but there is a function WSModelManipulator,I have tried it, it also get a error.

Another small bug in modem.Data()

The second argument name in the call to mtz.Z in function read_data_file is incorrect, near line 1125 of modem.py:

z_err_array=z_dummy.copy().real
should be
zerr_array=z_dummy.copy().real

(same with tipper_err_array, should be tippererr_array)

p.s. and other references to z_errand err_tipperin the same function (should be zerrand errtipper)

Cheers
-Ralf

MT 2d anisotropic forward modeling

I would like to know if this tool can do the MT 2d anisotropic forward modeling and what are steps I can follow to produce the correct results?

Tabs vs spaces

I see Jared added a whole bunch of tab characters the other day to setup.py (not sure where else)

It still runs for me, but that's just lucky. All advice is not to mix spaces with tabs, and pick one or the other and stick with it. http://stackoverflow.com/questions/119562/tabs-versus-spaces-in-python-programming

What are we using for indentation?

Duplicate versions of MTPlotTools.py

There are two versions of MTPlotTools.py -- one in core and the other in imaging.

One needs to be deleted, but I'm not sure what code is duplicated over both, and whether other modules refer to one or other or both. Jared?

Bugs in Occam2D classes

I'm having quite a bit of trouble getting Occam2D classes to read in my existing Occam 2D runs.

Problems:

1. Windows paths with colons in the startup/iter files break the read2Diter method. I fixed that in at least the one place where I found it.
2. Code assumes model files are called INMODEL, and the keyword argument provided for it in the constructor for Occam2DModel basically does nothing. I fixed this for at least my case...
3. read2DInmodel() doesn't properly read the valid model file I'm using

The last one is the problem I haven't been able to fix. I don't know the algorithm used and haven't been able to sort it out yet. The model file is here: https://gist.github.com/3997973#file_model.txt

WinglinkTools only for 2D? Extension for 3D?

It took me a while now, where I wanted to test the WIngLinkTools and they just did not work. Is it correct, that they are specifically written for WingLink2D? I was confused by the comments 'coordinates X,Y,Z' at some point in the codes.

If so, what is needed for the 3D output then?

Can not run example script `birrp_2in2out_simple.py`

Can any of contributors demonstrate (possibly add a readme?) how to run the examples?
I tried this:

python birrp_2in2out_simple.py birrp/birrp_exec bp02/133 bp02/133

In header = MTfh.read_ts_header(fn) we have this line:
headerlements = ['samplingrate', 't_min', 'nsamples', 'lat', 'lon', 'elev']

None of these strings ['samplingrate', 't_min', 'nsamples', 'lat', 'lon', 'elev'] exist anywhere in the examples folder. So the headers for the data are not found at all, and the python script exits with this error:

No MTpy data files for station BP02/133 found in directory /home/sudipta/Dropbox/GA/mtpy/examples/data/EarthDataLogger/BP02/133 !!

Export Occam 2D xyz

Hi, Does anyone know of a way in the MT-py codes to export an ASCII file (with x,y,z,resistivity) from the outputs of an Occam inversion (rather than simply displaying it as an image)? Ideally it would be in real world coordinates.
Cheers,
(sorry is this is the wrong location to ask this question)

Occam2D graphic interface

HI,

I can't run the Occam2D graphic interface. Are there any modifications to do in the python's modules (especially the matplotlib module) ?
Thank you

Aboubakr

Python 2.4 Compatibility

I've had a bug report from a potential user that MTpy is not compatible with Python 2.4 or 2.5. It doesn't seem like a huge deal to make it 2.4 compatible (see error list below, mostly change with statement and logic syntax in lists) just some syntax changes. The questions is should we fork a branch that is compatible with 2.4 which will be compatible with 2.x or always use 2.4 syntax, which defeats some the great things about current Python versions. Any thoughts? I made a branch called Python_v_2.4 as a start to see how much work it would be. I'm also not positive about how compatible some of the major packages (Numpy, Scipy, Matplotlib, etc) we use are backwards compatible to 2.4.

running install
running build
running build_py
running build_scripts
running install_lib
byte-compiling /usr/lib/python2.4/site-packages/mtpy/utils/misc.py to misc.pyc
File "/usr/lib/python2.4/site-packages/mtpy/utils/misc.py", line 145
class MemoryCheck():
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/utils/configfile.py to configfile.pyc
File "/usr/lib/python2.4/site-packages/mtpy/utils/configfile.py", line 119
with open(filename) as F:
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/utils/remove_instrumentresponse_from_files.py to remove_instrumentresponse_from_files.pyc
File "/usr/lib/python2.4/site-packages/mtpy/utils/remove_instrumentresponse_from_files.py", line 81
lo_channels = list(set([i.upper() if len(i)==2 else 'B'+i.upper() for i in sys.argv[4].split(',')]))
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/utils/convert_coordinates_in_edis.py to convert_coordinates_in_edis.pyc
File "/usr/lib/python2.4/site-packages/mtpy/utils/convert_coordinates_in_edis.py", line 78
with open(infile,'r') as F:
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/utils/filehandling.py to filehandling.pyc
File "/usr/lib/python2.4/site-packages/mtpy/utils/filehandling.py", line 578
with open(tsfile,'r') as F:
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/test/utils.py to utils.pyc
File "/usr/lib/python2.4/site-packages/mtpy/test/utils.py", line 40
with open(i,'r') as f:
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/test/basic.py to basic.pyc
File "/usr/lib/python2.4/site-packages/mtpy/test/basic.py", line 23
with self.assertRaises(ValueError):
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/analysis/distortion.py to distortion.pyc
File "/usr/lib/python2.4/site-packages/mtpy/analysis/distortion.py", line 345
lo_dims = [ 4 if i == 1 else i for i in lo_dims ]
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/analysis/pt.py to pt.pyc
File "/usr/lib/python2.4/site-packages/mtpy/analysis/pt.py", line 922
class ResidualPhaseTensor():
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/processing/birrp.py to birrp.pyc
File "/usr/lib/python2.4/site-packages/mtpy/processing/birrp.py", line 239
with open(string_file,'w') as F:
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/processing/io.py to io.pyc
File "/usr/lib/python2.4/site-packages/mtpy/processing/io.py", line 58
with open(filename, 'r') as f:
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/modeling/occam1d.py to occam1d.pyc
File "/usr/lib/python2.4/site-packages/mtpy/modeling/occam1d.py", line 1803
class PlotL2():
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/modeling/occam2d.py to occam2d.pyc
File "/usr/lib/python2.4/site-packages/mtpy/modeling/occam2d.py", line 73
class Setup():
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/core/edi.py to edi.pyc
File "/usr/lib/python2.4/site-packages/mtpy/core/edi.py", line 176
with open(infile,'r') as F:
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/usgs/zonge.py to zonge.pyc
File "/usr/lib/python2.4/site-packages/mtpy/usgs/zonge.py", line 34
class ZongeMTFT():
^
SyntaxError: invalid syntax
byte-compiling /usr/lib/python2.4/site-packages/mtpy/usgs/zen.py to zen.pyc
File "/usr/lib/python2.4/site-packages/mtpy/usgs/zen.py", line 2137
class ZenBIRRP():
^
SyntaxError: invalid syntax
running install_scripts
changing mode of /usr/bin/ws2vtk.py to 755

Apparent resistivity and phase sections and maps

Hi

This is Roshan here I am unable to plot the apparent resistivity and phase sections and maps.
When I am trying to plot the same for the examples data files given in the manual it doesn't produce a similar image.
Please help.

Roshan

Numpy -> VTK conversion

Stephan mentioned, there is some conversion from Numpy into VTK for the paraview visualisation around? Is it the wsmt_pv (thus specifically for th ws output)?

A generalised version would be a good thing to have in 'utils'. I start trying to browse through it today and see, if I get an idea.

Bug for phase tensor ellipse

Hi all,

(very quiet here, is MTPy still developed?)

there seems to be a bug when plotting phase tensor ellipses using the skew angle (beta) with the 'skew_seg' colours.

BUG: negative angles don't get correctly mapped into the assigned color segments

Here is the relevant code in sub get_plot_color in module mtcolors.py:

for bb in range(bounds.shape[0]):                     
    if colorx >= bounds[bb] and colorx < bounds[bb+1]:
        cvar = float(bounds[bb])/bounds.max()         
        return get_color(cvar, cmap)                  

for negative value bounds[bb] is the upper bracket, which however should be the lower bracket.

Proposed solution - use positive values and multiply with the original sign (not thoroughly tested!)

        # Get the sign of the skew angle:                                               
        skew_sign = np.sign(colorx)                                                     
        # New variable to work with                                                     
        abs_color = np.abs(colorx)                                                      
                                                                                        
        # Use original code, start with angles outside given range:                     
        if abs_color > bounds.max():                                                    
            # Restore sign here:                                                        
            cvar = skew_sign*1.0                                                        
            return get_color(cvar, cmap)                                                
                                                                                        
        # Take care of angles in interval:
        for bb in range(bounds.shape[0]):                                               
            if abs_color >= bounds[bb] and abs_color < bounds[bb+1]:                    
                # Restore sign here:                                                    
                cvar = skew_sign*float(bounds[bb])/bounds.max()                         
                return get_color(cvar, cmap)                                            

Cheers
-Ralf

New Occam2D module

I have just updated the Occam2D module with the new class structure. However, much functionality is still missing compared to the old occamtools (moved to legacy directory).
Any issues with the old occam functionalities, or can I just try to mimic them in the newer style?
Esp. the masking of points....??

Converion of EDI from Spectra to Impedance

Hi,
Is there any script to convert an EDI file in format of Spectra to Impedance?
Thanks,
Sasha

ModEM module

Dear all,
I can not run any scripts relate to ModEM and this error happens:
python ModEM_plot_slices.py
ERROR:root:the config yaml file /usr/local/lib/python2.7/dist-packages/mtpy-1.0.dev0-py2.7.egg/mtpy/utils/logging.yml does not exist?
None
Traceback (most recent call last):
File "ModEM_plot_slices.py", line 6, in
from mtpy.modeling.modem import PlotSlices
ImportError: No module named modem

If I want to make some changes in default code, is it possible ?

If yes, then how can I do that?
I want to make the y-axis in the phase section plot between the range of (-180, 180) with the spacing of 30.
I changed the major_locator(MultipleLocator) from 15 to 30 and saved that file.
But after running the program I did not get the outcome according to the changes.

setup.py

Just trying to run the setup.py and got a warning install_lib: 'build\lib' does not exist --no python modules to install. And nothing installed. was something changed?

Documentation

As we go we should start documenting what already works. We can wait till we're further along, or not, but it would be good to have some overviews, or examples.

I'm familiar with Sphinx, and know how to set it up so that we can have a website hosted by Github Pages, for example like http://emolch.github.com/pyrocko/. The advantage of doing it this way is that we can automatically retrieve docstrings into the HTML docs, to get pages like http://emolch.github.com/pyrocko/v0.3/trace.html

Or, it could be simpler -- we could just put some code examples on the wiki that's part of this repository already.

Most projects I've seen go with Sphinx, it's more flexible.

What do you think?

Python | Mtpy compatibilty

Which version of python is more likely to work?
I used the 2.7 but i had syntax errors (parentheses etc...)

bug in _set_zrot()

print 'length of angle list not correct'+
'({0} instead of {1})'.format(len(angle), len(self.z))

should be

print 'length of angle list not correct'+
'({0} instead of {1})'.format(len(angle), len(self.Z.z))

Rotation in legacy code

Hi all,

I encountered the following code in function _get_station_data (Class Occam2DData in legacy):

self.thetar = thetar
if abs(self.thetar) > 2*np.pi:
    self.thetar = self.thetar*(np.pi/180)

this converts degree-angles larger than 2PI (why is that?) to radian-angles

However all subsequent code seem to expect degree-angles.

Here is code just before the rotation is done (in function rotatematrix_incl_errors in utils/calculator.py):

phi = math.radians(degreeangle)

which again converts to radian.

I don't think this is correct.

Cheers
-Ralf

Small bug in modem

Hi all,

near line 699 in modem_new.py:

self.data_array[ii]['z_err'][jj] = interp_z.z_err[kk, :, :]

but interp_z has attribute zerrnot z_err

(same with interp_t.tippererr)

Cheers
-Ralf

Bug occam2d_rewrite

Hi all,

I encountered the following in the occam2d_rewrite module, near line 1058:

self.geoelectric_strike = np.median(np.nonzero(strike_angles))

this is a bug: np.nonzero return an array of integers (indices where strike_angles is non-zero).

Question: is this the correct module to use for creating an Occam2D datafile?

Cheers
-Ralf

ImportError: No module named LatLongUTMconversion

C:\Users\nietky\github\MTpy [master]> python -c 'import MTpy.core.MTtools'
Traceback (most recent call last):
  File "<string>", line 1, in <module>
  File "MTpy\core\MTtools.py", line 14, in <module>
    import Z
  File "MTpy\core\Z.py", line 18, in <module>
    import LatLongUTMconversion as utm2ll
ImportError: No module named LatLongUTMconversion

Make sure you use the full package path on imports.

New MT Py at https://github.com/MTgeophysics/mtpy

MTPy is now being further developed by Geoscience Australia, at

https://github.com/MTgeophysics/mtpy

Thanks GA!