Commit 5a3db32f authored by Jim King's avatar Jim King Committed by Jørgen Lind
Browse files

poststack passes tests, prestack offset and trace flag isn't correct

parent 5547b8b3
......@@ -45,7 +45,13 @@ bool
SEGYFileInfo::HasGatherOffset() const
{
// TODO unbinned gathers may be angle gathers?
return m_segyType == SEGY::SEGYType::Prestack || IsUnbinned();
return m_segyType == SEGYType::Prestack || m_segyType == SEGYType::Prestack2D || IsUnbinned();
}
bool
SEGYFileInfo::Is2D() const
{
return m_segyType == SEGYType::Poststack2D || m_segyType == SEGYType::Prestack2D;
}
SEGYBinInfo
......
......@@ -117,6 +117,8 @@ struct SEGYFileInfo
OPENVDS_EXPORT bool IsUnbinned() const;
OPENVDS_EXPORT bool HasGatherOffset() const;
OPENVDS_EXPORT bool Is2D() const;
};
#endif // SEGY_FILE_INFO_H
......@@ -60,8 +60,8 @@ def test_2d_poststack_volume_info(poststack_2d_segy, output_vds):
assert descriptor.coordinateMax == 1977.0
assert descriptor.coordinateStep == 1.0
assert layout.getChannelValueRangeMin(0) == -302681.94
assert layout.getChannelValueRangeMax(0) == 303745.38
assert layout.getChannelValueRangeMin(0) == -302681.9375
assert layout.getChannelValueRangeMax(0) == 303745.375
# check lattice metadata is NOT there
assert not layout.isMetadataDoubleVector2Available(
......@@ -103,6 +103,21 @@ def test_2d_poststack_read(poststack_2d_segy, output_vds):
assert result == 0, ex.output()
assert Path(output_vds.filename).exists()
with openvds.open(output_vds.filename, "") as handle:
layout = openvds.getLayout(handle)
access_manager = openvds.getAccessManager(handle)
dim0_size = layout.getDimensionNumSamples(0)
dim1_size = layout.getDimensionNumSamples(1)
request = access_manager.requestVolumeSubset((0, 0, 0, 0, 0, 0), (dim0_size, dim1_size, 1, 1, 1, 1),
dimensionsND=openvds.DimensionsND.Dimensions_01,
format=openvds.VolumeDataChannelDescriptor.Format.Format_R32)
data = request.data.reshape(dim1_size, dim0_size)
for dim1 in range(dim1_size):
total = 0
for dim0 in range(dim0_size):
total += abs(data[dim1, dim0])
assert total > 0.0, f"trace at {dim1}"
def test_2d_prestack_volume_info(prestack_2d_segy, output_vds):
# run importer
......@@ -126,26 +141,6 @@ def test_2d_prestack_volume_info(prestack_2d_segy, output_vds):
assert result == 0, ex.output()
assert Path(output_vds.filename).exists()
def test_2d_prestack_read(prestack_2d_segy, output_vds):
# run importer
# open VDS
# read VDS (entire thing at once)
# check that each trace is non-zero
ex = ImportExecutor()
ex.add_arg("--2d")
ex.add_arg("--prestack")
ex.add_args(["--vdsfile", output_vds.filename])
ex.add_arg(prestack_2d_segy)
result = ex.run()
assert result == 0, ex.output()
assert Path(output_vds.filename).exists()
with openvds.open(output_vds.filename, "") as handle:
layout = openvds.getLayout(handle)
......@@ -177,8 +172,8 @@ def test_2d_prestack_read(prestack_2d_segy, output_vds):
assert descriptor.coordinateMax == 1977.0
assert descriptor.coordinateStep == 1.0
assert layout.getChannelValueRangeMin(0) == -420799.12
assert layout.getChannelValueRangeMax(0) == 421478.38
assert layout.getChannelValueRangeMin(0) == -420799.125
assert layout.getChannelValueRangeMax(0) == 421478.375
# check lattice metadata is NOT there
assert not layout.isMetadataDoubleVector2Available(
......@@ -200,3 +195,69 @@ def test_2d_prestack_read(prestack_2d_segy, output_vds):
openvds.KnownMetadata.energySourcePointNumbers().name)
assert layout.isMetadataBLOBAvailable(openvds.KnownMetadata.ensembleNumbers().category,
openvds.KnownMetadata.ensembleNumbers().name)
def test_2d_prestack_read(prestack_2d_segy, output_vds):
# run importer
# open VDS
# read VDS (entire thing at once)
# check that each trace is non-zero
ex = ImportExecutor()
ex.add_arg("--2d")
ex.add_arg("--prestack")
ex.add_args(["--vdsfile", output_vds.filename])
ex.add_arg(prestack_2d_segy)
result = ex.run()
assert result == 0, ex.output()
assert Path(output_vds.filename).exists()
with openvds.open(output_vds.filename, "") as handle:
layout = openvds.getLayout(handle)
access_manager = openvds.getAccessManager(handle)
dim0_size = layout.getDimensionNumSamples(0)
dim1_size = layout.getDimensionNumSamples(1)
dim2_size = layout.getDimensionNumSamples(2)
trace_channel = layout.getChannelIndex("Trace")
assert trace_channel > 0
offset_channel = layout.getChannelIndex("Offset")
assert offset_channel > 0
request = access_manager.requestVolumeSubset((0, 0, 0, 0, 0, 0),
(dim0_size, dim1_size, dim2_size, 1, 1, 1),
format=openvds.VolumeDataChannelDescriptor.Format.Format_R32,
dimensionsND=openvds.DimensionsND.Dimensions_012)
offset_request = access_manager.requestVolumeSubset((0, 0, 0, 0, 0, 0),
(1, dim1_size, dim2_size, 1, 1, 1),
format=openvds.VolumeDataChannelDescriptor.Format.Format_R32,
dimensionsND=openvds.DimensionsND.Dimensions_012)
trace_flag_request = access_manager.requestVolumeSubset((0, 0, 0, 0, 0, 0),
(1, dim1_size, dim2_size, 1, 1, 1),
format=openvds.VolumeDataChannelDescriptor.Format.Format_U8,
dimensionsND=openvds.DimensionsND.Dimensions_012)
data = request.data.reshape(dim2_size, dim1_size, dim0_size)
offset_data = offset_request.data.reshape(dim2_size, dim1_size)
trace_flag_data = trace_flag_request.data.reshape(dim2_size, dim1_size)
for dim2 in range(dim2_size):
for dim1 in range(dim1_size):
trace_flag = trace_flag_data[dim2, dim1]
trace_offset = offset_data[dim2, dim1]
ensemble_number = layout.getAxisDescriptor(2).sampleIndexToCoordinate(dim2)
total = 0
for dim0 in range(dim0_size):
total += abs(data[dim1, dim0])
message = f"trace at dim1 {dim1} dim2 {dim2} offset {trace_offset} CDP {ensemble_number}"
if trace_flag == 0:
assert total == 0.0, f"dead {message}"
else:
assert total > 0.0, message
add_executable(SEGYImport
SEGYImport.cpp
TraceInfo2DManager.cpp
TraceInfo2DManager.h
)
target_link_libraries(SEGYImport PRIVATE Threads::Threads openvds segyutils jsoncpp_lib_static fmt::fmt)
......
This diff is collapsed.
/****************************************************************************
** Copyright 2021 The Open Group
** Copyright 2021 Bluware, Inc.
**
** Licensed under the Apache License, Version 2.0 (the "License");
** you may not use this file except in compliance with the License.
** You may obtain a copy of the License at
**
** http://www.apache.org/licenses/LICENSE-2.0
**
** Unless required by applicable law or agreed to in writing, software
** distributed under the License is distributed on an "AS IS" BASIS,
** WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
** See the License for the specific language governing permissions and
** limitations under the License.
****************************************************************************/
#include "TraceInfo2DManager.h"
TraceInfo2DManagerStore::TraceInfo2DManagerStore(const TraceInfo2D* pData, int64_t count) : m_traceInfo(pData, pData + count)
{
}
void TraceInfo2DManagerStore::Clear()
{
m_traceInfo.clear();
ClearEnsembleNumberMap();
}
int TraceInfo2DManagerStore::Count()
{
return static_cast<int>(m_traceInfo.size());
}
const TraceInfo2D& TraceInfo2DManagerStore::Get(int traceNumber)
{
return m_traceInfo.at(traceNumber);
}
const TraceInfo2D& TraceInfo2DManagerStore::Front()
{
return m_traceInfo.front();
}
int TraceInfo2DManagerStore::GetIndexOfEnsembleNumber(int ensembleNumber)
{
if (m_ensembleToIndex.size() != m_traceInfo.size())
{
BuildEnsembleNumberMap();
}
const auto
iter = m_ensembleToIndex.find(ensembleNumber);
assert(iter != m_ensembleToIndex.end());
if (iter == m_ensembleToIndex.end())
{
return 0;
}
return iter->second;
}
void TraceInfo2DManagerStore::BuildEnsembleNumberMap()
{
m_ensembleToIndex.clear();
for (int index = 0; index < m_traceInfo.size(); ++index)
{
m_ensembleToIndex.insert(std::make_pair(m_traceInfo[index].ensembleNumber, index));
}
}
void TraceInfo2DManagerStore::ClearEnsembleNumberMap()
{
m_ensembleToIndex.clear();
}
TraceInfo2DManagerImpl::TraceInfo2DManagerImpl(SEGY::Endianness endianness, double scaleOverride, const SEGY::HeaderField& coordinateScaleHeaderField, const SEGY::HeaderField& xCoordinateHeaderField, const SEGY::HeaderField& yCoordinateHeaderField, const SEGY::HeaderField& startTimeHeaderField, const SEGY::HeaderField& espnHeaderField, const SEGY::HeaderField& ensembleNumberHeaderField)
: m_endianness(endianness), m_scaleOverride(scaleOverride), m_coordinateScaleHeaderField(coordinateScaleHeaderField), m_xCoordinateHeaderField(xCoordinateHeaderField), m_yCoordinateHeaderField(yCoordinateHeaderField), m_startTimeHeaderField(startTimeHeaderField), m_espnHeaderField(espnHeaderField), m_ensembleNumberHeaderField(ensembleNumberHeaderField), m_previousEnsembleNumber(-1)
{
}
void TraceInfo2DManagerImpl::Clear()
{
TraceInfo2DManagerStore::Clear();
m_previousEnsembleNumber = -1;
}
void TraceInfo2DManagerImpl::AddTraceInfo(const char* traceHeader)
{
double coordScaleFactor = 1.0;
if (m_scaleOverride == 0.0)
{
int scale = SEGY::ReadFieldFromHeader(traceHeader, m_coordinateScaleHeaderField, m_endianness);
if (scale < 0)
{
coordScaleFactor = 1.0 / static_cast<double>(std::abs(scale));
}
else if (scale > 0)
{
coordScaleFactor = 1.0 * static_cast<double>(std::abs(scale));
}
}
else
{
coordScaleFactor = m_scaleOverride;
}
TraceInfo2D traceInfo{};
traceInfo.energySourcePointNumber = SEGY::ReadFieldFromHeader(traceHeader, m_espnHeaderField, m_endianness);
traceInfo.ensembleNumber = SEGY::ReadFieldFromHeader(traceHeader, m_ensembleNumberHeaderField, m_endianness);
traceInfo.x = coordScaleFactor * SEGY::ReadFieldFromHeader(traceHeader, m_xCoordinateHeaderField, m_endianness);
traceInfo.y = coordScaleFactor * SEGY::ReadFieldFromHeader(traceHeader, m_yCoordinateHeaderField, m_endianness);
traceInfo.startTime = SEGY::ReadFieldFromHeader(traceHeader, m_startTimeHeaderField, m_endianness);
if (m_previousEnsembleNumber != traceInfo.ensembleNumber)
{
m_traceInfo.push_back(traceInfo);
m_previousEnsembleNumber = traceInfo.ensembleNumber;
// any modification to trace info collection should invalidate the ensemble number mapping
ClearEnsembleNumberMap();
}
}
/****************************************************************************
** Copyright 2021 The Open Group
** Copyright 2021 Bluware, Inc.
**
** Licensed under the Apache License, Version 2.0 (the "License");
** you may not use this file except in compliance with the License.
** You may obtain a copy of the License at
**
** http://www.apache.org/licenses/LICENSE-2.0
**
** Unless required by applicable law or agreed to in writing, software
** distributed under the License is distributed on an "AS IS" BASIS,
** WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
** See the License for the specific language governing permissions and
** limitations under the License.
****************************************************************************/
#pragma once
#include "SEGYUtils/SEGYFileInfo.h"
#include <unordered_map>
struct TraceInfo2D
{
int energySourcePointNumber;
int ensembleNumber;
double x;
double y;
int startTime;
};
class TraceInfo2DManager
{
public:
virtual ~TraceInfo2DManager() = default;
virtual void Clear() {}
virtual void AddTraceInfo(const char* traceHeader) {}
virtual int Count() { return 0; }
virtual const TraceInfo2D& Get(int traceNumber) { return m_placeHolder; }
virtual const TraceInfo2D& Front() { return m_placeHolder; }
virtual int GetIndexOfEnsembleNumber(int ensembleNumber) { return 0; }
private:
TraceInfo2D m_placeHolder{};
};
class TraceInfo2DManagerStore : public TraceInfo2DManager
{
public:
TraceInfo2DManagerStore() = default;
TraceInfo2DManagerStore(const TraceInfo2D* pData, int64_t count);
void Clear() override;
int Count() override;
const TraceInfo2D& Get(int traceNumber) override;
const TraceInfo2D& Front() override;
int GetIndexOfEnsembleNumber(int ensembleNumber) override;
protected:
std::vector<TraceInfo2D> m_traceInfo;
std::unordered_map<int, int>
m_ensembleToIndex;
void BuildEnsembleNumberMap();
void ClearEnsembleNumberMap();
};
class TraceInfo2DManagerImpl : public TraceInfo2DManagerStore
{
public:
TraceInfo2DManagerImpl(SEGY::Endianness endianness, double scaleOverride, const SEGY::HeaderField& coordinateScaleHeaderField, const SEGY::HeaderField& xCoordinateHeaderField, const SEGY::HeaderField& yCoordinateHeaderField, const SEGY::HeaderField& startTimeHeaderField, const SEGY::HeaderField& espnHeaderField, const SEGY::HeaderField& ensembleNumberHeaderField);
void Clear() override;
void AddTraceInfo(const char* traceHeader) override;
private:
SEGY::Endianness m_endianness;
double m_scaleOverride;
SEGY::HeaderField m_coordinateScaleHeaderField;
SEGY::HeaderField m_xCoordinateHeaderField;
SEGY::HeaderField m_yCoordinateHeaderField;
SEGY::HeaderField m_startTimeHeaderField;
SEGY::HeaderField m_espnHeaderField;
SEGY::HeaderField m_ensembleNumberHeaderField;
int m_previousEnsembleNumber;
};
\ No newline at end of file
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment