Errors to save 3D numpy array in segy file

[victortxa]

I'm trying to write in a segy file a 3D data that I calculated and it is in the format of numpy.ndarray with shape (501, 64860) (first dimension is time), so each column is a trace.
First I tried to follow the example Make segy file from sctrach (https://github.com/Statoil/segyio#examples):

import segyio
import numpy as np
spec = segyio.spec()
filename = '/home/victor/output/output_file.sgy'

spec.sorting = 2
spec.format = 1
spec.samples = 501 #time sample number
spec.ilines = np.arange(345) #inline number
spec.xlines = np.arange(188) #xline number

with segyio.create(filename, spec) as f:

    # write the line itself to the file and the inline number in all this line's headers
    for ilno in spec.ilines:
        f.iline[ilno] = f.iline[ilno] = np.zeros(
            (188, 501), dtype=np.single) + ilno
        f.iline[ilno] = d5_out.T[ilno*len(spec.xlines):len(spec.xlines)+ilno*len(spec.xlines),:]
        f.header.iline[ilno] = {
            segyio.TraceField.INLINE_3D: ilno,
            segyio.TraceField.offset: 0
        }

    # then do the same for xlines
    for xlno in spec.xlines:
        f.header.xline[xlno] = {
           segyio.TraceField.CROSSLINE_3D: xlno,
           segyio.TraceField.TRACE_SAMPLE_INTERVAL: 2000
        }

Then, it raised the error:

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-136-a7139c546f3c> in <module>()
----> 1 with segyio.create(filename, spec) as f:
      2 
      3     # write the line itself to the file and the inline number in all this line's headers
      4     for ilno in spec.ilines:
      5         f.iline[ilno] = f.iline[ilno] = np.zeros(

~/anaconda3/lib/python3.6/site-packages/segyio/create.py in create(filename, spec)
    144     _segyio.putfield(binary, 3213, tracecount)
    145     _segyio.putfield(binary, 3217, 4000)
--> 146     _segyio.putfield(binary, 3221, len(samples))
    147     _segyio.putfield(binary, 3225, int(spec.format))
    148     _segyio.putfield(binary, 3505, int(ext_headers))

TypeError: len() of unsized object

My second attempt was following (in some way) this:

spec = segyio.spec()
filename = '/home/victor/output/test4.sgy'

spec.sorting = 2
spec.format = 1
spec.samples = list(range(501))
spec.ilines = list(range(345))
spec.xlines = list(range(188))
spec.offsets = range(188)
with segyio.create(filename, spec) as f:
    nx, no = len(spec.xlines), len(spec.offsets)
    
    # write the line itself to the file and the inline number in all this line's headers
    for ilno in spec.ilines:
        for xlno in spec.xlines:
            ix = xlno + ilno*len(spec.xlines)
            f.trace[ix] = d5_out.T[ix, :]
            f.header[ix] = { segyio.TraceField.INLINE_3D: ilno,
                             segyio.TraceField.CROSSLINE_3D: xlno,
                             segyio.TraceField.offset: xlno}

I know there is at least one problem in this second attempt, the variable spec.offsets, the reason for the error:

---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
<ipython-input-161-18a82df13c3f> in <module>()
     16             ix = xlno + ilno*len(spec.xlines) #(ilindex * nx * no) + (tr * no) + offset_index
---> 17             f.trace[ix] = d5_out.T[ix,:]
     18             f.header[ix] = { segyio.TraceField.INLINE_3D: ilno,
     19                              segyio.TraceField.CROSSLINE_3D: xlno,

~/anaconda3/lib/python3.6/site-packages/segyio/_trace.py in __setitem__(self, index, val)
     46 
     47         if not 0 <= abs(index) < len(self):
---> 48             raise IndexError(self.index_errmsg.format(index, len(self)-1))
     49 
     50         self.write_trace(index, val, self._file)

IndexError: Trace 3984 not in range [-3983,3983]

If I don't use this variable (as if I don't have the info to write in the header, so that I comment #spec.offsets = range(188)), the error message is:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-165-2dfbb450c12e> in <module>()
     16             ix = xlno + ilno*len(spec.xlines) #(ilindex * nx * no) + (tr * no) + offset_index
---> 17             f.trace[ix] = d5_out.T[ix,:]
     18             f.header[ix] = { segyio.TraceField.INLINE_3D: ilno,
     19                              segyio.TraceField.CROSSLINE_3D: xlno}#,

~/anaconda3/lib/python3.6/site-packages/segyio/_trace.py in __setitem__(self, index, val)
     45             return
     46 
---> 47         if not 0 <= abs(index) < len(self):
     48             raise IndexError(self.index_errmsg.format(index, len(self)-1))
     49 

ValueError: __len__() should return >= 0

Just to mention, when I run the example from the repository main page (Make segy file from sctrach), a similar error message is shown:

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-166-91270423ed28> in <module>()
     11 spec.xlines = np.arange(20)
     12 
---> 13 with segyio.create(filename, spec) as f:
     14 
     15     # write the line itself to the file and the inline number in all this line's headers

~/anaconda3/lib/python3.6/site-packages/segyio/create.py in create(filename, spec)
    144     _segyio.putfield(binary, 3213, tracecount)
    145     _segyio.putfield(binary, 3217, 4000)
--> 146     _segyio.putfield(binary, 3221, len(samples))
    147     _segyio.putfield(binary, 3225, int(spec.format))
    148     _segyio.putfield(binary, 3505, int(ext_headers))

TypeError: len() of unsized object

[jokva]

Hi, and thanks for the detailed report!

First of all, I'm sorry, but the example in the https://github.com/Statoil/segyio#examples section is slightly outdated, since samples isn't a list. I'll have that updated. You seem to have figured that out yourself, though, which is good. offsets needs to be a list too - are you sure you have 188 of them?

Presumably, your volume is post-stack (in which offsets should have a length of 1), and your spec should be

spec.samples = list(range(501))
spec.ilines = list(range(345))
spec.xlines = [1]

edit: corrected it should be:

spec.samples = list(range(501))
spec.ilines = list(range(345))
spec.xlines = list(range(188))
spec.offsets = [1]

[victortxa]

For each inline I have 188 xlines (each trace of a inline is part of a xline), because of this I wrote spec.xlines = list(range(188)). Shouldn't I do this? This is the same reason I wrote spec.offsets = range(188).
If I wrote spec.xlines = [1] I think it would be as the data is 2D since only the inline number is changing. In other words:

ix=1 , xlno=1 , ilno=0
ix=2 , xlno=1 , ilno=1
ix=3 , xlno=1 , ilno=2
ix=4 , xlno=1 , ilno=3
ix=5 , xlno=1 , ilno=4
ix=6 , xlno=1 , ilno=5
ix=7 , xlno=1 , ilno=6
...

IndexError: Trace 345 not in range [-344,344]

The range shown above is the double of the length of all inlines (actually 2*ilno-1 in the range [-344,344]).

I think I'm making some confusion with the values...however the error seems related to some range values that I don't know from where it comes (Trace 3984 not in range [-3983,3983]). Do you have any clues?

[jokva]

spec.xlines = list(range(188)) is correct and perfectly fine. offsets however, are only present in pre-stack cubes (and your inline-crossline ranges are consistent with post-stack).

If you do spec.xlines = [1] the data would indeed be 2D.

345*188 gives 64860 and matches your ndarray perfectly. I realise now that my code was completely 100% wrong because I meant to also copy your xline range, and only replace your offsets. Try this:

spec.sorting = 2
spec.format = 1
spec.samples = range(501) #time sample number
spec.ilines = np.arange(345) #inline number
spec.xlines = np.arange(188) #xline number
spec.offsets = [1]

spec.offsets needs to be set to a range, usually of 1 element. This may be defaulted in the future.

[victortxa]

Unfortunately it raises another error:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-185-36ddfe2331fc> in <module>()
     16             ix = xlno + ilno*len(spec.xlines) #(ilindex * nx * no) + (tr * no) + offset_index
     17             print('ix={} , xlno={} , ilno={}'.format(ix,xlno,ilno))
---> 18             f.trace[ix] = d5_out.T[ix,:]
     19             f.header[ix] = { segyio.TraceField.INLINE_3D: ilno,
     20                              segyio.TraceField.CROSSLINE_3D: xlno}#,

~/anaconda3/lib/python3.6/site-packages/segyio/_trace.py in __setitem__(self, index, val)
     45             return
     46 
---> 47         if not 0 <= abs(index) < len(self):
     48             raise IndexError(self.index_errmsg.format(index, len(self)-1))
     49 

ValueError: __len__() should return >= 0

It still seem related to the dimensions, but I can't see the problem!

[jokva]

Interesting. That one should show up if the tracecount is not set up properly.

What does print(f.tracecount) give you?

[jokva]

For what it's worth, I just ran this:

>>> import segyio
>>> spec = segyio.spec()
>>> spec.sorting = 2
>>> spec.format = 1
>>> spec.samples = list(range(501))
>>> spec.ilines = np.arange(345)
>>> spec.xlines = np.arange(188)
>>> spec.offsets = [1]
>>> f = segyio.create('file.sgy', spec)
>>> f.trace[345] = np.zeros(501)

No errors.

[victortxa]

It is really weird: it raised the same error and the output of tracecount doesn't make sense to me (-676). Shouldn't it be 501?

Could be something related to my python installation or maybe another library I use?

-676

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-5-025a327f6bd6> in <module>()
      8 f = segyio.create('/home/victor/output/file.sgy', spec)
      9 print(f.tracecount)
---> 10 f.trace[0] = np.zeros(501)

~/anaconda3/lib/python3.6/site-packages/segyio/_trace.py in __setitem__(self, index, val)
     45             return
     46 
---> 47         if not 0 <= abs(index) < len(self):
     48             raise IndexError(self.index_errmsg.format(index, len(self)-1))
     49 

ValueError: __len__() should return >= 0

[jokva]

I doubt there's something installed, but you never know. What's your segyio version? Can you post your full create function and all spec arguments?

This does give me a few ideas for usability improvements and error messages.

[victortxa]

The version is '1.5.1'.
It follows the complete function:

import segyio
import numpy as np

spec = segyio.spec()
filename = '/home/victor/output/d5.sgy'

spec.sorting = 2
spec.format = 1
spec.samples = list(range(501))
spec.ilines = np.arange(345)
spec.xlines = np.arange(188)
spec.offsets = [1]

#****************FIRST attempt*************
with segyio.create(filename, spec) as f:

    # write the line itself to the file and the inline number in all this line's headers
    for ilno in spec.ilines:
        f.iline[ilno] = f.iline[ilno] = np.zeros(
            (188, 501), dtype=np.single) + ilno
        f.iline[ilno] = d5_out.T[ilno*188:188+ilno*188,:]
        f.header.iline[ilno] = {
            segyio.TraceField.INLINE_3D: ilno,
            segyio.TraceField.offset: 0
        }

    # then do the same for xlines
    for xlno in spec.xlines:
        f.header.xline[xlno] = {
           segyio.TraceField.CROSSLINE_3D: xlno,
           segyio.TraceField.TRACE_SAMPLE_INTERVAL: 2000
        }

---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
<ipython-input-13-a7139c546f3c> in <module>()
      8         f.header.iline[ilno] = {
      9             segyio.TraceField.INLINE_3D: ilno,
---> 10             segyio.TraceField.offset: 0
     11         }
     12 

~/anaconda3/lib/python3.6/site-packages/segyio/_line.py in __setitem__(self, lineno, val)
    139         else:
    140             t0 = self._index(lineno, offset)
--> 141             self.writefn(t0, self.len, self.stride, val)
    142 
    143     def __len__(self):

~/anaconda3/lib/python3.6/site-packages/segyio/_header.py in writefn(self, t0, length, stride, val)
     74 
     75         for i, x in zip(range(start, stop, stride), val):
---> 76             self[i] = x
     77 
     78     @property

~/anaconda3/lib/python3.6/site-packages/segyio/_header.py in __setitem__(self, traceno, val)
     44             traceno = traceno[0]
     45 
---> 46         self.__getitem__(traceno, buf).update(val)
     47 
     48     def __iter__(self):

~/anaconda3/lib/python3.6/site-packages/segyio/_header.py in __getitem__(self, traceno, buf)
     33             return gen()
     34 
---> 35         return Field.trace(buf, traceno=traceno, segy=self.segy)
     36 
     37     def __setitem__(self, traceno, val):

~/anaconda3/lib/python3.6/site-packages/segyio/_field.py in trace(cls, buf, traceno, segy)
    100 
    101         if traceno >= segy.tracecount or traceno < 0:
--> 102             raise IndexError("Header out of range: 0 <= {} < {}".format(traceno, segy.tracecount))
    103 
    104         if buf is None:

IndexError: Header out of range: 0 <= 0 < -676

#****************SECOND attempt*************
spec = segyio.spec()
filename = '/home/victor/output/test4.sgy'

spec.sorting = 2
spec.format = 1
spec.samples = list(range(501))
spec.ilines = np.arange(345)
spec.xlines = np.arange(188)#list(range(188))
spec.offsets = [1]#range(188)
with segyio.create(filename, spec) as f:
    nx, no = len(spec.xlines), len(spec.offsets)
    
    # write the line itself to the file and the inline number in all this line's headers
    for ilno in spec.ilines:
        for xlno in spec.xlines:
            ix = xlno + ilno*len(spec.xlines) #(ilindex * nx * no) + (tr * no) + offset_index
            f.trace[ix] = d5_out.T[ix,:]
            f.header[ix] = { segyio.TraceField.INLINE_3D: ilno,
                             segyio.TraceField.CROSSLINE_3D: xlno,
                             segyio.TraceField.offset: xlno}
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-17-597c95265356> in <module>()
     15         for xlno in spec.xlines:
     16             ix = xlno + ilno*len(spec.xlines) #(ilindex * nx * no) + (tr * no) + offset_index
---> 17             f.trace[ix] = d5_out.T[ix,:]
     18             f.header[ix] = { segyio.TraceField.INLINE_3D: ilno,
     19                              segyio.TraceField.CROSSLINE_3D: xlno,

~/anaconda3/lib/python3.6/site-packages/segyio/_trace.py in __setitem__(self, index, val)
     45             return
     46 
---> 47         if not 0 <= abs(index) < len(self):
     48             raise IndexError(self.index_errmsg.format(index, len(self)-1))
     49 

ValueError: __len__() should return >= 0
            

I did a quick test and maybe I found the problem...I think it is the amount of traces/file size. I cropped my data and then a file was generated without error message (using the "second method" to write the file). I was slightly increasing the size but it seems there is a limit of 42,190,800 bytes in the file size (when I tried to read a error message said the size indicated in the header was different form the actual data size).
I'll verify the integrity of the file and make others tests.

[jokva]

Interesting, what system are you running this on. This really looks like a 32/64 bit issue from your description, and I'd love to get to the bottom of that.

[victortxa]

I'm using Ubuntu 4.13.0-36-generic #40~16.04.1-Ubuntu and Python 3.6.3 :: Anaconda custom (64-bit). I don't know the details of segyio implementation, but it seems size of data is also a problem in obspy for writing segy files.

[jokva]

Hi,

I just reproduced your error on my Debian Stretch python 3.5 implementation, so this is not a system-specific issue for you. I'll figure out a fix for you!

[victortxa]

Thank you!!
Just other 2 comments I forgot:
I just tested to read the file also using obspy, but it complains. I don't know how reliable obspy can be, but for the file generated from segyio it says: Too little data left in the file to unpack it according to its trace header. This is most likely either due to a wrong byte order or a corrupt file.
The other observation is: I tested increasing the number of inline recorded to verify the limit for the segy size. Although I could generate the file until about 31960 traces, when reading the file (with segyio) it showed that only about 18000 traces was recorded, and the rest was zeros. Maybe I did somethinh wrong, but it can be a clue about the data size limitation.

[jokva]

Ok, I think I figured it out. It was an overflow, but due to a 16-bit integer.

Some context: when I read the standard, I must've misinterpreted (and really glossed over a VERY important detail) the wording "Number of data traces per ensemble" to mean "traces-in-this-file". This field is 2 bytes and overflows right away.

During a refactoring pass, in order to generalise, I used this field to pass size information from create down to the primitives, and you just hit this limit. The segyio version I tested against must've been too old or something, so I couldn't reproduce right away. Your code was perfectly fine all along (save samples being an array, but we figured that out). For this, I thank you for your patience.

I'll issue a fix on this as soon as possible.

[victortxa]

Very good news!!
Thank you for your time and patience too!! .segy always give headaches for everyone.

[jokva]

Ok, so I've written a tiny patch that should fix this issue. #236

I'm gonna let the CI do its business before I merge it into master (it should be quite safe) and release the next bugfix release of segyio. If you're savvy enough to build segyio yourself from source, I'd love it if you checked out the changes and tried if it solves your problem. If it does, the bugfix should be freely available sometime tomorrow.

[jokva]

I have to run now, but I just want to thank you for using segyio and many thanks for a great bug report! We're obviously happy to help you with any other issues you might have.

[victortxa]

I've cloned and built it successfully, but I don't remember how to 'force' importing the version I've build.
I'll try to do.

[jokva]

The typical way to do that would be either to overwrite your existing install (assuming it's pip --user it would be in ~/.local, otherwise typically in /usr/local/) or manipulate the PYTHONPATH environment variable to make sure your new version is picked up.

[victortxa]

My problem is Python is installed via Anaconda, so only editing the ` PYTHONPATH` didn't work... I'll try a work around asap.

…

On Wed, Mar 14, 2018 at 8:59 PM, Jørgen Kvalsvik ***@***.***> wrote: The typical way to do that would be either to overwrite your existing install (assuming it's pip --user it would be in ~/.local, otherwise typically in /usr/local/) or manipulate the PYTHONPATH environment variable to make sure your new version is picked up. — You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub <#235 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AJrEtjbFJjVi68GdN2rfVF8jj7E5F12Xks5tea74gaJpZM4SpRuy> .

-- Att, Victor Almeida.

[victortxa]

With my user and Anaconda installation I didn't find how to 'force' the build version.
So I did again in another machine, however it seems something was not successful. There wasn't any error message during sudo make install.

---------------------------------------------------------------------------
ImportError                               Traceback (most recent call last)
<ipython-input-1-3b14e2f7415e> in <module>()
----> 1 import segyio

/usr/local/lib/python2.7/dist-packages/segyio/__init__.pyc in <module>()
     86 from .su import su
     87 from .open import open
---> 88 from .create import create
     89 from .segy import SegyFile, spec
     90 from .tools import dt, sample_indexes, create_text_header, native

/usr/local/lib/python2.7/dist-packages/segyio/create.py in <module>()
      2 import numpy
      3 import segyio
----> 4 import segyio._segyio as _segyio
      5 
      6 

ImportError: libsegyio.so.1: cannot open shared object file: No such file or directory

[jokva]

Is this on Ubuntu? The installation script might be broken, and not properly installing the library to /usr/local/lib64, in which case you won't find it. If you have a look in /usr/local/lib and /usr/local/lib64, that might be a clue.

[victortxa]

Yes, the file is in /usr/local/lib

[jokva]

Ok, that means that the cmake GNUInstallDirs doesn't work properly.

You can manually copy or symlink the file from /usr/local/lib/ to /usr/local/lib64 to make it work. I'll look into the cmake thing.

[victortxa]

Manually copying didn't work, but creating symlink in lib64 to lib worked!

[victortxa]

I just tested now and it seems the error is still there:

In [1]: import segyio

In [2]: import numpy as np

In [3]: d5_out=np.zeros([501,188*345],dtype=np.uint32)

In [4]: spec = segyio.spec()
   ...: filename = '/home/lagex/test4.sgy'

In [5]: spec.sorting = 2
   ...: spec.format = 1
   ...: spec.samples = list(range(501))
   ...: spec.ilines = np.arange(345)
   ...: spec.xlines = np.arange(188)#list(range(188))
   ...: spec.offsets = [1]#range(188)
   ...: 

In [6]: with segyio.create(filename, spec) as f:
   ...:     nx, no = len(spec.xlines), len(spec.offsets)
   ...:     
   ...:     # write the line itself to the file and the inline number in all this line's headers
   ...:     for ilno in spec.ilines:
   ...:         for xlno in spec.xlines:
   ...:             ix = xlno + ilno*len(spec.xlines) #(ilindex * nx * no) + (tr * no) + offset_index
   ...:             f.trace[ix] = d5_out.T[ix,:]
   ...:             f.header[ix] = { segyio.TraceField.INLINE_3D: ilno,
   ...:                              segyio.TraceField.CROSSLINE_3D: xlno,
   ...:                              segyio.TraceField.offset: xlno}
   ...:             
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-6-11138f6f85e6> in <module>()
      6         for xlno in spec.xlines:
      7             ix = xlno + ilno*len(spec.xlines) #(ilindex * nx * no) + (tr * no) + offset_index
----> 8             f.trace[ix] = d5_out.T[ix,:]
      9             f.header[ix] = { segyio.TraceField.INLINE_3D: ilno,
     10                              segyio.TraceField.CROSSLINE_3D: xlno,

/usr/local/lib/python2.7/dist-packages/segyio/_trace.pyc in __setitem__(self, index, val)
     45             return
     46 
---> 47         if not 0 <= abs(index) < len(self):
     48             raise IndexError(self.index_errmsg.format(index, len(self)-1))
     49 

ValueError: __len__() should return >= 0

Is it related to the trace length? I remember in somewhere it complained about a negative value for some input/output it was using.

[jokva]

No, it's a bug internal to create - it triggers a 16-bit integer overflow.

Did you apply my patch, or is this a clean checkout?

[victortxa]

I built from your repository (https://github.com/jokva/segyio). Is there with the patch, right?

[jokva]

Only if you check-out the create-tracecount-arg-not-in-header branch. :--)

[victortxa]

hehe...
Sorry about that! I just read the "Fix bug" in the master branch...Testing now!

[victortxa]

@jokva I tested and the file was written! I write and read with segyio without problems!
I also tried to read with obspy but it complained about some inconsistency between header and data. However it's possible I didn't write all the correct parameters for obspy read properly.
Later I'll test to read in OpendTect. I let you know how it's going.

[jokva]

Not surprisng that obspy breaks - by default, only the absolute minimum for segyio to read it is written, and that even excludes the trace headers, which must be written manually.

Great that it now works - this has surely been a useful report. I'll figure out how to test for this and merge the patch upstream. Thank you again for your contribution.

Please let me know how it works with OpendTect.

[sumitav]

I was asking how to convert a numpy array(as data) to segy file.Previously have converted the same to segy file but it's 2d in nature i.e., I can't create 3d parameters say depth at particular inline/crossline,,can't define inline and cross line of my own etc etc..
So given a 3d numpy array how to convert that to a seg-y file with necessary headers and data.???
Can anyone please solve this

[jokva]

https://segyio.readthedocs.io/en/latest/segyio.html#segyio.tools.from_array3D should do what you want.

[sumitav]

@jokva have seen it but the error I am now facing tells ":AttributeError: module 'segyio' has no attribute 'tools'" on running this command"segyio.tools.from_array3D(path,data)".
I did tried to solve it but couldn't.'
Can you please solve this please @victortxa @asbjorn @ingvald @jokva @markusdregi @lightenup @epa095 ..

[sumitav]

@jokva got the seg-y file from 3d array but it doesn't contain any data ,,like while plotting it I can't access the specified inline and crossline also when am viewing it in my segy tool its void.
Kindly help me out