Merge pull request #413 from jokva/unsigned-16bit-header-word

Interpret some header words as unsigned

applications/CMakeLists.txt

@@ -75,6 +75,7 @@ if (NOT BUILD_TESTING)
 endif ()
 
 set(small ${testdata}/small.sgy)
+set(long  ${testdata}/long.sgy)
 set(test ${CMAKE_CURRENT_SOURCE_DIR}/test)
 add_test(NAME catr.arg.t1       COMMAND segyio-catr -t 0 ${small})
 add_test(NAME catr.arg.t2       COMMAND segyio-catr -t 5 ${small})
@@ -122,6 +123,7 @@ add_custom_command(
     COMMAND segyio-catb ${small} > catb.out
     COMMAND segyio-catb ${small} -n -d > catbnd.out
     COMMAND segyio-cath ${small} > cath.out
+    COMMAND segyio-catb ${long} > catb-long.out
     COMMAND segyio-catr -r 4 9 -t 12 ${small} > catr.out
     COMMAND segyio-catr -t 1 ${small} -n -d > catrnd.out
 )
@@ -141,6 +143,10 @@ add_test(NAME catb.output.nd
          COMMAND ${CMAKE_COMMAND} -E compare_files ${test}/catbnd.output
                                                    catbnd.out
 )
+add_test(NAME catb-long.output
+         COMMAND ${CMAKE_COMMAND} -E compare_files ${test}/catb-long.output
+                                                   catb-long.out
+)
 add_test(NAME cath.output
          COMMAND ${CMAKE_COMMAND} -E compare_files ${test}/cath.output cath.out
 )

applications/segyio-catb.c

@@ -26,6 +26,17 @@ static int printhelp(){
 static int get_binary_value( char* binheader, int bfield ){
     int32_t f;
     segy_get_bfield( binheader, bfield, &f );
+
+    /*
+     * convert cannot-be-negative values to unsigned int, as mandated by SEGY-Y
+     * rev2
+     */
+    switch (bfield) {
+        case SEGY_BIN_SAMPLES:
+        case SEGY_BIN_SAMPLES_ORIG:
+            f = (int32_t)((uint16_t)(f));
+            break;
+    }
     return f;
 }
 

applications/segyio-catr.c

@@ -696,6 +696,17 @@ int main( int argc, char** argv ) {
                 int f;
                 segy_get_field( trheader, fields[j], &f );
                 if( opts.nonzero && !f ) continue;
+
+                /*
+                 * convert cannot-be-negative values to unsigned int, as mandated by SEGY-Y
+                 * rev2
+                 */
+                switch (f) {
+                    case SEGY_TR_SAMPLE_COUNT:
+                        f = (int)((uint16_t)f);
+                        break;
+                }
+
                 if( opts.description ) {
                     printf( "%s\t%d\t%d\t%s\n",
                             labels[j],

applications/test/catb-long.output

@@ -0,0 +1,30 @@
+jobid	0
+lino	0
+reno	0
+ntrpr	1
+nart	1
+hdt	1000
+dto	1000
+hns	60000
+nso	43392
+format	1
+fold	0
+tsort	0
+vscode	0
+hsfs	0
+hsfe	0
+hslen	0
+hstyp	0
+schn	0
+hstas	0
+hstae	0
+htatyp	0
+hcorr	0
+bgrcv	0
+rcvm	0
+mfeet	0
+polyt	0
+vpol	0
+rev	0
+trflag	0
+exth	0

lib/CMakeLists.txt

@@ -70,6 +70,7 @@ if (NOT BUILD_TESTING)
 endif ()
 
 configure_file(${testdata}/small.sgy                          test-data/small.sgy                            COPYONLY)
+configure_file(${testdata}/long.sgy                           test-data/long.sgy                             COPYONLY)
 configure_file(${testdata}/small-ps-dec-il-xl-off.sgy         test-data/small-ps-dec-il-xl-off.sgy           COPYONLY)
 configure_file(${testdata}/small-ps-dec-il-inc-xl-off.sgy     test-data/small-ps-dec-il-inc-xl-off.sgy       COPYONLY)
 configure_file(${testdata}/small-ps-dec-xl-inc-il-off.sgy     test-data/small-ps-dec-xl-inc-il-off.sgy       COPYONLY)

lib/include/segyio/segy.h

@@ -49,7 +49,9 @@ int segy_binheader( segy_file*, char* buf );
 int segy_write_binheader( segy_file*, const char* buf );
 /*
  * exception: the int returned is the number of samples (the segy standard only
- * allocates 2 octets for this, so it comfortably sits inside an int
+ * allocates 2 octets for this, so it comfortably sits inside an int. If the
+ * value is negative, it is converted to an unsigned 16-bit integer, allowing
+ * up to 65536 samples per trace.
  */
 int segy_samples( const char* binheader );
 /*

lib/src/segy.c

@@ -968,6 +968,7 @@ int segy_set_format( segy_file* fp, int format ) {
 int segy_samples( const char* binheader ) {
     int32_t samples = 0;
     segy_get_bfield( binheader, SEGY_BIN_SAMPLES, &samples );
+    samples = (int32_t)((uint16_t)samples);
     return samples;
 }
 

lib/test/segy.cpp

@@ -1355,6 +1355,28 @@ SCENARIO( "reading text header", "[c.segy]" ) {
         CHECK( ascii == expected );
 }
 
+TEST_CASE("open file with >32k traces", "[c.segy]") {
+    const char* file = "test-data/long.sgy";
+    unique_segy ufp(segy_open(file, "rb"));
+    auto fp = ufp.get();
+
+    char bin[SEGY_BINARY_HEADER_SIZE];
+    Err err = segy_binheader(fp, bin);
+    REQUIRE(err == Err::ok());
+
+    /*
+     * This file was created with 60000 samples and three traces. Before SEG-Y
+     * rev2 this had a max-size of 2^15 - 1 because it was 2s complement
+     * 16-bit. However, rev2 allows this to be unsigned 16-bit as a negative
+     * sample-size is non-sensical.
+     *
+     * segy_samples is aware of this, but the get_bfield function returns
+     * signed/unsigned faithfully
+     */
+    const auto samples = segy_samples(bin);
+    CHECK(samples == 60000);
+}
+
 SCENARIO( "reading a large file", "[c.segy]" ) {
     GIVEN( "a large file" ) {
         const char* file = "4G-file.sgy";

python/segyio/segyio.cpp

@@ -396,6 +396,9 @@ PyObject* segycreate( segyiofd* self, PyObject* args, PyObject* kwargs ) {
     if( samples <= 0 )
         return ValueError( "expected samples > 0" );
 
+    if (samples > std::numeric_limits< std::uint16_t >::max())
+        return ValueError("samples must fit in an uint16");
+
     if( tracecount <= 0 )
         return ValueError( "expected tracecount > 0" );
 
@@ -1440,6 +1443,18 @@ PyObject* getfield( PyObject*, PyObject *args ) {
             : segy_get_field(  buffer.buf< const char >(), field, &value )
             ;
 
+    /*
+     * Some fields be negative, and SEG-Y rev2 considers them unsigned. When
+     * reading them, widen to unsigned.
+     */
+    switch (field) {
+        case SEGY_TR_SAMPLE_COUNT:
+        case SEGY_BIN_SAMPLES:
+        case SEGY_BIN_SAMPLES_ORIG:
+            value = int(std::uint16_t(value));
+            break;
+    }
+
     switch( err ) {
         case SEGY_OK:            return PyLong_FromLong( value );
         case SEGY_INVALID_FIELD: return KeyError( "No such field %d", field );

python/test/segy.py

@@ -1483,6 +1483,33 @@ def test_create_sgy_skip_lines(tmpdir):
         assert 3.014 == approx(f.iline[3][1][4], abs=1e-4)
         assert 7.023 == approx(f.iline[7][2][3], abs=1e-4)
 
+def test_create_overflowing_samples_fail(tmpdir):
+    spec = segyio.spec()
+    spec.tracecount = 1
+    spec.samples = list(range(np.power(2, 16)))
+    spec.format = 1
+    with pytest.raises(ValueError):
+        with segyio.create(tmpdir / '65k.sgy', spec) as _:
+            pass
+
+def test_create_uint16_samples(tmpdir):
+    spec = segyio.spec()
+    spec.tracecount = 1
+    spec.samples = list(range(np.power(2, 16) - 1))
+    spec.format = 1
+
+    with segyio.create(tmpdir / '65k.sgy', spec) as f:
+        assert len(f.samples) == 65535
+        assert f.bin[segyio.su.hns] == len(f.samples)
+        f.trace[0] = spec.samples
+        f.header[0] = {
+            segyio.su.ns: len(f.samples),
+        }
+
+    with segyio.open(tmpdir / '65k.sgy', ignore_geometry = True) as f:
+        assert len(f.samples) == 65535
+        assert f.bin[segyio.su.hns] == len(f.samples)
+        assert f.header[0][segyio.su.ns] == len(f.samples)
 
 def mklines(fname):
     spec = segyio.spec()