1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
|
From: =?utf-8?q?Picca_Fr=C3=A9d=C3=A9ric-Emmanuel?=
<picca@synchrotron-soleil.fr>
Date: Tue, 4 May 2021 10:35:03 +0200
Subject: fix for mrcimage reader
---
fabio/mrcimage.py | 19 +++++++++----
fabio/test/codecs/test_mrcimage.py | 57 +++++++++++++++++++++-----------------
2 files changed, 45 insertions(+), 31 deletions(-)
diff --git a/fabio/mrcimage.py b/fabio/mrcimage.py
index 2cb95b1..f379aee 100644
--- a/fabio/mrcimage.py
+++ b/fabio/mrcimage.py
@@ -31,13 +31,16 @@ email: Jerome.Kieffer@terre-adelie.org
Specifications from:
http://ami.scripps.edu/software/mrctools/mrc_specification.php
+New version on:
+https://www.ccpem.ac.uk/mrc_format/mrc_format.php
+https://www.fei-software-center.com/tem-apps/MRC-2014-Specifications/
"""
__authors__ = ["Jérôme Kieffer"]
__contact__ = "Jerome.Kieffer@terre-adelie.org"
__license__ = "MIT"
__copyright__ = "Jérôme Kieffer"
-__date__ = "03/04/2020"
+__date__ = "23/04/2021"
import logging
import numpy
@@ -85,12 +88,18 @@ class MrcImage(FabioImage):
int_block = numpy.frombuffer(infile.read(56 * 4), dtype=numpy.int32)
for key, value in zip(self.KEYS, int_block):
self.header[key] = value
- if self.header["MAP"] != 542130509:
- logger.info("Expected 'MAP ', got %s", self.header["MAP"].tobytes())
+ # convert some headers ...
+ self.header["MAP"] = self.header["MAP"].tobytes().decode()
+ if self.header["MAP"][:3] not in ('MAP ', 'FEI'):
+ logger.info("Expected 'MAP ', got %s", self.header["MAP"])
for i in range(10):
label = "LABEL_%02i" % i
- self.header[label] = infile.read(80).strip()
+ self.header[label] = infile.read(80).decode().strip()
+
+ # Read extended header
+ self.header["extended"] = infile.read(self.header["NSYMBT"])
+
dim1 = int(self.header["NX"])
dim2 = int(self.header["NY"])
self._shape = dim2, dim1
@@ -126,7 +135,7 @@ class MrcImage(FabioImage):
:param frame: frame number
"""
assert frame < self.nframes
- return 1024 + frame * self.imagesize
+ return 1024 + self.header["NSYMBT"] + frame * self.imagesize
def _makeframename(self):
self.filename = "%s$%04d" % (self.sequencefilename,
diff --git a/fabio/test/codecs/test_mrcimage.py b/fabio/test/codecs/test_mrcimage.py
index c2a84a6..5643045 100755
--- a/fabio/test/codecs/test_mrcimage.py
+++ b/fabio/test/codecs/test_mrcimage.py
@@ -46,43 +46,48 @@ from ..utilstest import UtilsTest
class TestMrc(unittest.TestCase):
"""basic test"""
- mrcfilename = 'EMD-3001.map'
- npyfilename = 'EMD-3001.npy'
- results = """EMD-3001.map 73 43 -0.36814222 0.678705 0.006804995 0.1630334"""
+ mrcfilename = ('EMD-3001.map', "0p67_5s_0000.mrc")
+ npyfilename = ('EMD-3001.npy', "0p67_5s_0000.npy")
+ results = """EMD-3001.map 73 43 -0.36814222 0.678705 0.0062340125 0.16349247
+ 0p67_5s_0000.mrc 2048 2048 -344.0 33553.0 82.653259 243.213061"""
def setUp(self):
"""Download files"""
self.fn = {}
- for i in [self.mrcfilename, self.npyfilename]:
+ for i in self.mrcfilename + self.npyfilename:
self.fn[i] = UtilsTest.getimage(i + ".bz2")[:-4]
for i in self.fn:
assert os.path.exists(self.fn[i])
def test_read(self):
- """ check we can read kcd images"""
- vals = self.results.split()
- dim1, dim2 = [int(x) for x in vals[1:3]]
- shape = dim2, dim1
- mini, maxi, mean, stddev = [float(x) for x in vals[3:]]
- for ext in ["", ".gz", ".bz2"]:
- try:
- obj = openimage(self.fn[self.mrcfilename] + ext)
- except Exception as err:
- logger.error("unable to read: %s", self.fn[self.mrcfilename] + ext)
- raise err
- self.assertAlmostEqual(mini, obj.getmin(), 4, "getmin" + ext)
- self.assertAlmostEqual(maxi, obj.getmax(), 4, "getmax" + ext)
- self.assertAlmostEqual(mean, obj.getmean(), 4, "getmean" + ext)
- self.assertAlmostEqual(stddev, obj.getstddev(), 4, "getstddev" + ext)
- self.assertEqual(shape, obj.shape, "shape" + ext)
+ """ check we can read mrc images"""
+ for results in self.results.split("\n"):
+ vals = results.split()
+ dim1, dim2 = [int(x) for x in vals[1:3]]
+ shape = dim2, dim1
+ mini, maxi, mean, stddev = [float(x) for x in vals[3:]]
+ for ext in ["", ".gz", ".bz2"]:
+ filename = self.fn[vals[0]] + ext
+ try:
+ obj = openimage(filename)
+ except Exception as err:
+ logger.error("unable to read: %s", filename)
+ raise err
+ self.assertAlmostEqual(mini, obj.getmin(), 4, f"{filename} getmin")
+ self.assertAlmostEqual(maxi, obj.getmax(), 4, f"{filename} getmax")
+ self.assertAlmostEqual(mean, obj.getmean(), 4, f"{filename} getmean")
+ self.assertAlmostEqual(stddev, obj.getstddev(), 4, f"{filename} getstddev")
+ self.assertEqual(shape, obj.shape, f"{filename} shape")
def test_same(self):
- """ see if we can read kcd images and if they are the same as the EDF """
- mrc = MrcImage()
- mrc.read(self.fn[self.mrcfilename])
- npy = fabio.open(self.fn[self.npyfilename])
- diff = (mrc.data - npy.data)
- self.assertAlmostEqual(abs(diff).sum(), 0, 4)
+ """ see if we can read mrc images and if they are the same as the numpy dump"""
+ for mrcfilename, npyfilename in zip(self.mrcfilename, self.npyfilename):
+ logger.info("Comparing files %s and %s", mrcfilename, npyfilename)
+ mrc = MrcImage()
+ mrc.read(self.fn[mrcfilename])
+ npy = fabio.open(self.fn[npyfilename])
+ diff = (mrc.data - npy.data)
+ self.assertAlmostEqual(abs(diff).sum(), 0, 4, msg=mrcfilename)
def suite():
|
