Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add read function for friction.data files and a test #557

Closed
wants to merge 1 commit into from
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
46 changes: 44 additions & 2 deletions src/python/geoclaw/data.py
Original file line number Diff line number Diff line change
Expand Up @@ -371,8 +371,6 @@ def write(self,data_source='setrun.py', out_file='dtopo.data'):
def read(self, path, force=False):
r"""Read a dtopography data file."""

print(self.dtopofiles)

with open(os.path.abspath(path), 'r') as data_file:

file_name = None
Expand Down Expand Up @@ -611,6 +609,7 @@ def __init__(self):
# File support
self.add_attribute('friction_files', [])


def write(self, out_file='friction.data', data_source='setrun.py'):

self.open_data_file(out_file, data_source)
Expand Down Expand Up @@ -647,6 +646,49 @@ def write(self, out_file='friction.data', data_source='setrun.py'):
self.close_data_file()


def read(self, path="friction.data"):

super(FrictionData, self).read(path)

with open(os.path.abspath(path), 'r') as data_file:
# Skip header
for line in data_file:
if "#" != line[0]:
break
value, tail = data_file.readline().split("=:")
self.variable_friction = bool(value)
value, tail = data_file.readline().split("=:")
self.friction_index = int(value)
data_file.readline()

value, tail = data_file.readline().split("=:")
num_friction_regions = int(value)
data_file.readline()

# Read in each region
self.friction_regions = []
for i in range(num_friction_regions):
# Read in region
region = []
for j in range(4):
values, tail = data_file.readline().split("=:")
region.append( [float(value) for value in values.split()] )
# Check
if len(region[3]) != len(region[2]) - 1:
raise ValueError("Incorrect number of depths and coefficients.")

self.friction_regions.append(region)
data_file.readline()

# Read in files
value, tail = data_file.readline().split("=:")
num_friction_files = int(value)
if num_friction_files > 0:
raise NotImplementedError("Friction files are not supported ",
"yet.")



class MultilayerData(clawpack.clawutil.data.ClawData):
r"""
Multilayer SWE data object
Expand Down
63 changes: 63 additions & 0 deletions tests/test_data.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,63 @@
#!/usr/bin/env python
# encoding: utf-8

"""Tests for reading and writing GeoClaw data files"""

import tempfile
import shutil
import os

import numpy

import clawpack.clawutil.test as test
import clawpack.geoclaw.data


def test_read_friction_data():
r"""Test readinga and writing of FrictionData files"""

# Test Data
test_regions = []
test_regions.append([(-99.0, -70.0), (8.0, 32.0),
[numpy.infty, 0.0, -numpy.infty],
[0.030, 0.022]])
test_regions.append([(-98, 25.25), (-90, 30),
[numpy.infty, -10.0, -200.0, -numpy.infty],
[0.030, 0.012, 0.022]])

# Create temp directory
temp_path = tempfile.mkdtemp()

try:
data_file = os.path.join(temp_path, "friction.data")

# Create and write out data object
friction_data = clawpack.geoclaw.data.FrictionData()
friction_data.variable_friction = True
for i in range(2):
friction_data.friction_regions.append(test_regions[i])
friction_data.write(data_file)

# Read data object
read_friction_data = clawpack.geoclaw.data.FrictionData()
read_friction_data.read(data_file)

# Tests
for (i, region) in enumerate(test_regions):
for j in range(4):
assert numpy.allclose(region[j],
read_friction_data.friction_regions[i][j])

except AssertionError as e:
# If the assertion failed then copy the contents of the directory
shutil.copytree(temp_path, os.path.join(os.getcwd(),
"test_read_friction_data"))
raise e

finally:
shutil.rmtree(temp_path)


if __name__ == "__main__":
test_read_friction_data()
print("All tests passed.")