MTgeophysics / mtpy

Python toolbox for standard Magnetotelluric (MT) data analysis
GNU General Public License v3.0
145 stars 66 forks source link

MT_obj not convertible to .edi file with only Tipper information #110

Open k-a-mendoza opened 4 years ago

k-a-mendoza commented 4 years ago

I made an Mt_obj with just tipper information:


from mtpy.core.mt import MT
def create_simple_metadata():
    """
    extracts metadata from the pandas dataframe row, specifically:
    site (lat, lon, elevation, station, line)
    acquisition (azimuth, recording time)
    """
    lat     = 0
    lon     = 0
    elev    = 0
    station ='1'

    return lat, lon, elev, station

def generate_tipper_arrays():
    """
    function to populate numpy arrays from pandas dataframe
    """
    frequencies = np.asarray([11, 22, 45, 90, 180, 360, 720, 1440])
    t_n=[]
    t_e=[]
    for freq in frequencies:
        t_e.append(0.5 + 0.5j)
        t_n.append(0.5 + 0.5j)

    t_n = np.asarray(t_n,dtype=np.complex64)
    t_e = np.asarray(t_e,dtype=np.complex64)
    return frequencies, t_n, t_e

def create_tipper():
    """
    creates an MtPy Tipper array from source numpy arrays and relevant metadata
    """
    frequencies, t_n, t_e = generate_tipper_arrays()

    array = np.zeros((len(frequencies),1,2),dtype=np.complex64)
    array[:,0,0]=t_n
    array[:,0,1]=t_e

    tipper = Tipper(tipper_array=array, tipper_err_array=None,
                    freq=frequencies)
    return tipper

def parse_row_into_mt_obj():
    lat, lon, elev, station = create_simple_metadata(row)
    tipper   = create_tipper(row)
    mt_obj   = MT(station=station,lat=lat,lon=lon,elev=elev,Tipper=tipper)
    return mt_obj

mt_obj = parse_row_into_mt_obj()
mt_obj.write_mt_file(save_dir='', fn_basename='1', file_type='edi',
                      longitude_format='LON',latlon_format='dd')

This is the error that pops up

/anaconda3/envs/mtpy_inv/lib/python3.6/site-packages/mtpy/core/mt.py in write_mt_file(self, save_dir, fn_basename, file_type, new_Z_obj, new_Tipper_obj, longitude_format, latlon_format)
    436                                       new_Tipper=new_Tipper_obj,
    437                                       longitude_format=longitude_format,
--> 438                                       latlon_format=latlon_format)
    439         elif file_type == 'xml':
    440             fn = self._write_xml_file(fn,

~/anaconda3/envs/mtpy_inv/lib/python3.6/site-packages/mtpy/core/mt.py in _write_edi_file(self, new_edi_fn, new_Z, new_Tipper, longitude_format, latlon_format)
    650 
    651         # get mtsec
--> 652         edi_obj.Data_sect = self._edi_set_data_sect()
    653 
    654         # set Z and T data

~/anaconda3/envs/mtpy_inv/lib/python3.6/site-packages/mtpy/core/mt.py in _edi_set_data_sect(self)
    885         sect = MTedi.DataSection()
    886         sect.data_type = 'MT'
--> 887         sect.nfreq = self.Z.z.shape[0]
    888         sect.sectid = self.station
    889         nchan = 5

AttributeError: 'NoneType' object has no attribute 'shape'

Is this expected behavior? Is there a recommended way to work with .edi files and live mt_obj's with only tipper data?

alkirkby commented 4 years ago

At the moment MTPy expects that there will be impedance information to write an edi file. There are a series of fixes that will need to be implemented to allow tipper only. We'll work on this. In the meantime, you could put a dummy z array in with null values to allow the code to work. Hope this helps.

kujaku11 commented 4 years ago

Agree with Alison, the simplest would be to make a dummy impedance tensor with the same number of frequencies as your induction vectors.