File: mrtrix_utils.h

package info (click to toggle)
mrtrix3 3.0.4-1
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid, trixie
  • size: 13,712 kB
  • sloc: cpp: 129,776; python: 9,494; sh: 593; makefile: 234; xml: 47
file content (192 lines) | stat: -rw-r--r-- 6,677 bytes parent folder | download
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
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
/* Copyright (c) 2008-2022 the MRtrix3 contributors.
 *
 * This Source Code Form is subject to the terms of the Mozilla Public
 * License, v. 2.0. If a copy of the MPL was not distributed with this
 * file, You can obtain one at http://mozilla.org/MPL/2.0/.
 *
 * Covered Software is provided under this License on an "as is"
 * basis, without warranty of any kind, either expressed, implied, or
 * statutory, including, without limitation, warranties that the
 * Covered Software is free of defects, merchantable, fit for a
 * particular purpose or non-infringing.
 * See the Mozilla Public License v. 2.0 for more details.
 *
 * For more details, see http://www.mrtrix.org/.
 */

#ifndef __formats_mrtrix_utils_h__
#define __formats_mrtrix_utils_h__

#include "header.h"
#include "types.h"

#include "file/gz.h"
#include "file/key_value.h"
#include "file/ofstream.h"
#include "file/path.h"
#include "file/utils.h"



namespace MR
{
  namespace Formats
  {


    // Read generic image header information - common between conventional, compressed and sparse formats
    template <class SourceType>
      void read_mrtrix_header (Header&, SourceType&);

    // These are helper functiosn for reading key/value pairs from either a File::KeyValue construct,
    //   or from a GZipped file (where the getline() function must be used explicitly)
    bool next_keyvalue (File::KeyValue::Reader&, std::string&, std::string&);
    bool next_keyvalue (File::GZ&,       std::string&, std::string&);

    // Get the path to a file - use same function for image data and sparse data
    // Note that the 'file' and 'sparse_file' fields are read in as entries in the map<string, string>
    //   by read_mrtrix_header(), and are therefore erased by this function so that they do not propagate
    //   into new images created using this header
    void get_mrtrix_file_path (Header&, const std::string&, std::string&, size_t&);

    // Write generic image header information to a stream -
    //   this could be an ofstream in the case of .mif, or a stringstream in the case of .mif.gz
    template <class StreamType>
      void write_mrtrix_header (const Header&, StreamType&);


    vector<ssize_t> parse_axes (size_t ndim, const std::string& specifier);




    template <class SourceType>
      void read_mrtrix_header (Header& H, SourceType& kv)
      {
        std::string dtype, layout;
        vector<uint64_t> dim;
        vector<default_type> vox, scaling;
        vector<vector<default_type>> transform;

        std::string key, value;
        while (next_keyvalue (kv, key, value)) {
          const std::string lkey = lowercase (key);
          if (lkey == "dim") dim = parse_ints<uint64_t> (value);
          else if (lkey == "vox") vox = parse_floats (value);
          else if (lkey == "layout") layout = value;
          else if (lkey == "datatype") dtype = value;
          else if (lkey == "scaling") scaling = parse_floats (value);
          else if (lkey == "transform")
            transform.push_back (parse_floats (value));
          else if (key.size() && value.size())
            add_line (H.keyval()[key], value); // Preserve capitalization if not a compulsory key
        }

        if (dim.empty())
          throw Exception ("missing \"dim\" specification for MRtrix image \"" + H.name() + "\"");
        H.ndim() = dim.size();
        for (size_t n = 0; n < dim.size(); n++) {
          if (dim[n] < 1)
            throw Exception ("invalid dimensions for MRtrix image \"" + H.name() + "\"");
          H.size(n) = dim[n];
        }

        if (vox.empty())
          throw Exception ("missing \"vox\" specification for MRtrix image \"" + H.name() + "\"");
        if (vox.size() < std::min (size_t(3), dim.size()))
          throw Exception ("too few entries in \"vox\" specification for MRtrix image \"" + H.name() + "\"");
        for (size_t n = 0; n < std::min<size_t> (vox.size(), H.ndim()); n++) {
          if (vox[n] < 0.0)
            throw Exception ("invalid voxel size for MRtrix image \"" + H.name() + "\"");
          H.spacing(n) = vox[n];
        }


        if (dtype.empty())
          throw Exception ("missing \"datatype\" specification for MRtrix image \"" + H.name() + "\"");
        H.datatype() = DataType::parse (dtype);


        if (layout.empty())
          throw Exception ("missing \"layout\" specification for MRtrix image \"" + H.name() + "\"");
        auto ax = parse_axes (H.ndim(), layout);
        for (size_t i = 0; i < ax.size(); ++i)
          H.stride(i) = ax[i];

        if (transform.size()) {

          auto check_transform = [&transform]() {
            if (transform.size() < 3) return false;
            for (auto row : transform)
              if (row.size() != 4)
                return false;
            return true;
          };
          if (!check_transform())
            throw Exception ("invalid \"transform\" specification for MRtrix image \"" + H.name() + "\"");

          for (int row = 0; row < 3; ++row)
            for (int col = 0; col < 4; ++col)
              H.transform() (row,col) = transform[row][col];
        }


        if (scaling.size()) {
          if (scaling.size() != 2)
            throw Exception ("invalid \"scaling\" specification for MRtrix image \"" + H.name() + "\"");
          H.set_intensity_scaling (scaling[1], scaling[0]);
        }

      }







    template <class StreamType>
      void write_mrtrix_header (const Header& H, StreamType& out)
      {
        out << "dim: " << H.size (0);
        for (size_t n = 1; n < H.ndim(); ++n)
          out << "," << H.size (n);

        out << "\nvox: " << H.spacing (0);
        for (size_t n = 1; n < H.ndim(); ++n)
          out << "," << H.spacing (n);

        auto stride = Stride::get (H);
        Stride::symbolise (stride);

        out << "\nlayout: " << (stride[0] >0 ? "+" : "-") << abs (stride[0])-1;
        for (size_t n = 1; n < H.ndim(); ++n)
          out << "," << (stride[n] >0 ? "+" : "-") << abs (stride[n])-1;

        DataType dt = H.datatype();
        dt.set_byte_order_native();
        out << "\ndatatype: " << dt.specifier();

        Eigen::IOFormat fmt(Eigen::FullPrecision, Eigen::DontAlignCols, ", ", "\ntransform: ", "", "", "\ntransform: ", "");
        out << H.transform().matrix().topLeftCorner(3,4).format(fmt);

        if (H.intensity_offset() != 0.0 || H.intensity_scale() != 1.0)
          out << "\nscaling: " << H.intensity_offset() << "," << H.intensity_scale();

        for (const auto& i : H.keyval())
          for (const auto& line : split_lines (i.second))
            out << "\n" << i.first << ": " << line;


        out << "\n";
      }





  }
}

#endif