« back to all changes in this revision
Viewing changes to pydendro/itrdb.py
- browse files at revision 35
- compare with another revision
- download diff
- download tarball
- view history from revision 35
- Committer: Brewster Malevich
- Date: 2012-03-14 17:05:30 UTC
- Revision ID: sbmalev@gmail.com-20120314170530-hed4r0qa1kb6tjj2
Name changes to prep for merge with other 'pydendro'.
- files added:
- pydendro
- files removed:
- dendro
added
removed
pydendro/itrdb.py
1
#! /usr/bin/env python
2
#
3
# Copyright 2011 S. Brewster Malevich <malevich@email.arizona.edu>
4
#
5
# PyDendro is free software: you can redistribute it and/or modify
6
# it under the terms of the GNU General Public License as published by
7
# the Free Software Foundation, either version 3 of the License, or
8
# (at your option) any later version.
9
#
10
# This program is distributed in the hope that it will be useful,
11
# but WITHOUT ANY WARRANTY; without even the implied warranty of
12
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13
# GNU General Public License for more details.
14
#
15
# You should have received a copy of the GNU General Public License
16
# along with this program. If not, see <http://www.gnu.org/licenses/>
17
18
"""
19
Classes to read dendrochronological files from the ITRDB.
20
"""
21
22
import os
23
import linecache
24
from scipy.stats.stats import pearsonr
25
import numpy as np
26
27
class DendroFile:
28
"""Parent class to be inherited specific dendro files formats.
29
30
Arguments:
31
path -- Path to the target .crn file.
32
33
"""
34
def __init__(self, path):
35
"""Opens file 'path' and checks that it exists.
36
"""
37
if not os.path.isfile(path):
38
raise IOError()
39
self._path = path
40
41
class RawMeas(DendroFile):
42
"""Class to handle *.rwl files from the ITRDB.
43
44
-This class is in progress.-
45
46
Arguments:
47
path -- Path to the target .crn file.
48
49
"""
50
def __init__(self, path):
51
"""Parses the RWL file, defining instance attributes.
52
"""
53
DendroFile.__init__(self,path)
54
# TODO: Test if file has header... (do they all?)
55
# TODO: Parse header file.
56
57
58
class Chronology(DendroFile):
59
"""Class to handle *.crn files from the ITRDB.
60
61
Arguments:
62
path -- Path to the target .crn file.
63
64
Returns:
65
An object of class Chronology.
66
67
68
Note that the optional ITRDB "Chronology Statistics" will not be
69
preserved. This object is based on format definitions of the ITRDB (see
70
http://www.ncdc.noaa.gov/paleo/treeinfo.html).
71
72
73
Instance attributes:
74
75
type
76
A string indicating the tree-ring chronology file type (e.g. "ARSTND",
77
"Standard", etc.).
78
79
years
80
A list of integers corresponding to each year with a chronology index
81
value.
82
83
values
84
A list of float for chronology tree-ring index values for each
85
of the chronology's years.
86
87
sampledepth
88
An integer list of the measurement sample site for each year of the
89
chronology.
90
91
speciesname
92
A string name of sampled species for the chronology site.
93
94
speciescode
95
A string code for the sampled species of the chronology site.
96
97
siteid
98
A string site ID for the chronology.
99
100
sitename
101
Name of the chronology's sample site. A string.
102
103
location
104
String giving the state and country of the chronology's sample site.
105
106
latlon
107
An integer of the chronology sample site's latitude and longitude
108
(ddmm or dddmm). I've found a variety of formats for this so, beware.
109
110
elevation
111
Integer of the chronology site's elevation.
112
113
firstyear
114
The first year of the chronology. An integer.
115
116
lastyear
117
The last year of the chronology. An integer.
118
119
collectiondate
120
A string giving the chonology's final collection date.
121
122
investigator
123
A string of lead investigator(s).
124
125
"""
126
def __init__(self, path):
127
"""Parses the Chronology file defining instance attributes.
128
"""
129
DendroFile.__init__(self,path)
130
self.type = None
131
self.years = []
132
self.values = []
133
self.sampledepth = []
134
self.speciesname = None
135
self.speciescode = None
136
self.siteid = None
137
self.sitename = None
138
self.location = None
139
self.latlon = None
140
self.elevation = None
141
self.firstyear = None
142
self.lastyear = None
143
self.collectiondate = None
144
self.investigator = None
145
victim = None
146
types = {'a.crn':'ARSTND',
147
'p.crn':'Low Pass Filter',
148
'r.crn':'Residual',
149
's.crn':'Standard',
150
'w.crn':'Re-Whitened Residual',
151
'n.crn':'Measurements Only'}
152
if self._path[-5:] in types:
153
self.type = types[self._path[-5:]]
154
else:
155
self.type = 'Unknown'
156
headermark = [' 1 ', ' 2 ', ' 3 '] # To identify the header.
157
header1 = linecache.getline(path, 1)
158
header2 = linecache.getline(path, 2)
159
header3 = linecache.getline(path, 3)
160
if header1[6:9] == headermark[0]:
161
self.siteid = str(header1[0:6]).strip()
162
self.sitename = str(header1[9:61]).strip()
163
self.speciescode = str(header1[61:65]).strip()
164
else:
165
raise NameError('File ' + self._path
166
+ ' appears to have no header.')
167
if header2[6:9] == headermark[1]:
168
self.location = str(header2[9:22]).strip()
169
self.speciesname = str(header2[22:30]).strip()
170
# Slice elevation to ':44' to cut out the "M".
171
self.elevation = int(header2[40:44])
172
self.latlon = str(header2[47:57]).strip()
173
self.firstyear = int(header2[67:71])
174
self.lastyear = int(header2[72:76])
175
self.years = range(self.firstyear, self.lastyear + 1)
176
else:
177
raise NameError('File ' + self._path
178
+ ' appears to have no header.')
179
if header3[6:9] == headermark[2]:
180
self.investigator = str(header3[9:72]).strip()
181
self.collectiondate = str(header3[72:80]).strip()
182
else:
183
raise NameError('File ' + self._path
184
+ ' appears to have no header.')
185
lastdecade = self.lastyear - (self.lastyear % 10) # To ID last line.
186
for line in open(path):
187
if line[6:9] in headermark:
188
continue
189
a = 10 # Tree-ring index index.
190
b = 14
191
c = 14 # Ring count index.
192
d = 17
193
for i in range(0, 10):
194
victim = line[a:b]
195
if victim == '9990':
196
a += 7
197
b += 7
198
c += 7
199
d += 7
200
continue
201
victim = victim.strip()
202
victim = victim.zfill(4)
203
self.values.append(float(victim[:1] + '.' + victim[1:]))
204
self.sampledepth.append(int(line[c:d]))
205
a += 7
206
b += 7
207
c += 7
208
d += 7
209
if int(line[6:10]) == lastdecade: # Stop after last line.
210
break
211
if len(self.values) != len(self.sampledepth) != len(self.years):
212
raise NameError('Unequal values, samples sizes'
213
+ 'and years from chronology file.')
214
215
def crop_corr(self, x, x_years):
216
"""Crop and correlate a chrnology with a given variable.
217
218
Arguments:
219
x -- An array/list(?) of your variable of interest.
220
x_years -- Array/list of years corresponding to the x observations.
221
222
Returns:
223
A tuple giving 1) Pearson's correlation coefficient, 2) 2-tailed
224
p-value, 3) length of the common period between crn and x.
225
226
"""
227
crn_values = np.array(self.values)
228
crn_years = np.array(self.years)
229
if len(crn_values) != len(crn_years):
230
raise ValueError('unequal values and years')
231
# Might want to check if x and x_years are not array...?
232
# Type checks are python taboo but I'm not sure how numpy will handle
233
# lists.
234
year_min = np.maximum(np.min(crn_years), np.min(x_years))
235
year_max = np.minimum(np.max(crn_years), np.max(x_years))
236
crn_crop = crn_values[(crn_years >= year_min)
237
& (crn_years <= year_max)]
238
x_crop = x[(x_years >= year_min) & (x_years <= year_max)]
239
results = pearsonr(crn_crop, x_crop)
240
results += (len(crn_crop),)
241
return results
242
Loggerhead is a web-based interface for Breezy
Version: 2.0.1