CiftiLib
A C++ library for CIFTI-2 and CIFTI-1 files
MatrixFunctions.h
1#ifndef __MATRIX_UTILITIES_H__
2#define __MATRIX_UTILITIES_H__
3
4/*LICENSE_START*/
5/*
6 * Copyright (c) 2014, Washington University School of Medicine
7 * All rights reserved.
8 *
9 * Redistribution and use in source and binary forms, with or without modification,
10 * are permitted provided that the following conditions are met:
11 *
12 * 1. Redistributions of source code must retain the above copyright notice,
13 * this list of conditions and the following disclaimer.
14 *
15 * 2. Redistributions in binary form must reproduce the above copyright notice,
16 * this list of conditions and the following disclaimer in the documentation
17 * and/or other materials provided with the distribution.
18 *
19 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
20 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO,
21 * THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
22 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
23 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
24 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
25 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
26 * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
27 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE,
28 * EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
29 */
30
31#include <vector>
32#include <cmath>
33#include "stdint.h"
34
35using namespace std;
36
37//NOTE: this is not intended to be used outside of FloatMatrix.cxx, error condition is a 0x0 matrix result
38//NOTE: if a matrix has a row shorter than the first row, expect problems. Calling checkDim will look for this, but it is not used internally, as FloatMatrix maintains this invariant
39
40namespace cifti {
41
43 {
44 typedef int64_t msize_t;//NOTE: must be signed due to using -1 as a sentinel
45
46 public:
50 template <typename T1, typename T2, typename T3, typename A>
51 static void multiply(const vector<vector<T1> > &left, const vector<vector<T2> > &right, vector<vector<T3> > &result);
52
56 template <typename T1, typename T2, typename T3>
57 static void multiply(const vector<vector<T1> > &left, const T2 right, vector<vector<T3> > &result);
58
62 template <typename T>
63 static void rref(vector<vector<T> > &inout);
64
68 template <typename T>
69 static void inverse(const vector<vector<T> > &in, vector<vector<T> > &result);
70
74 template <typename T1, typename T2, typename T3>
75 static void add(const vector<vector<T1> > &left, const vector<vector<T2> > &right, vector<vector<T3> > &result);
76
80 template <typename T1, typename T2, typename T3>
81 static void add(const vector<vector<T1> > &left, const T2 right, vector<vector<T3> > &result);
82
86 template <typename T1, typename T2, typename T3>
87 static void subtract(const vector<vector<T1> > &left, const vector<vector<T2> > &right, vector<vector<T3> > &result);
88
92 template <typename T>
93 static void transpose(const vector<vector<T> > &in, vector<vector<T> > &result);
94
98 template <typename T>
99 static bool checkDim(const vector<vector<T> > &in);
100
104 template <typename T>
105 static void resize(const msize_t rows, const msize_t columns, vector<vector<T> > &result, bool destructive = false);
106
110 template <typename T>
111 static void zeros(const msize_t rows, const msize_t columns, vector<vector<T> > &result);
112
116 template <typename T>
117 static void ones(const msize_t rows, const msize_t columns, vector<vector<T> > &result);
118
122 template <typename T>
123 static void identity(const msize_t size, vector<vector<T> > &result);
124
128 template <typename T1, typename T2, typename T3>
129 static void horizCat(const vector<vector<T1> > &left, const vector<vector<T2> > &right, vector<vector<T3> > &result);
130
134 template <typename T1, typename T2, typename T3>
135 static void vertCat(const vector<vector<T1> > &top, const vector<vector<T2> > &bottom, vector<vector<T3> > &result);
136
140 template <typename T>
141 static void getChunk(const msize_t firstrow, const msize_t lastrow, const msize_t firstcol, const msize_t lastcol, const vector<vector<T> > &in, vector<vector<T> > &result);
142
143 private:
147 template <typename T>
148 static void rref_big(vector<vector<T> > &inout);
149 };
150
151 template <typename T1, typename T2, typename T3, typename A>
152 void MatrixFunctions::multiply(const vector<vector<T1> >& left, const vector<vector<T2> >& right, vector<vector<T3> >& result)
153 {//the stupid multiply O(n^3) - the O(n^2.78) version might not be that hard to implement with the other functions here, but not as stable
154 msize_t leftrows = (msize_t)left.size(), rightrows = (msize_t)right.size(), leftcols, rightcols;
155 vector<vector<T3> > tempstorage, *tresult = &result;//pointer because you can't change a reference
156 bool copyout = false;
157 if (&left == &result || &right == &result)
158 {
159 copyout = true;
160 tresult = &tempstorage;
161 }
162 if (leftrows && rightrows)
163 {
164 leftcols = (msize_t)left[0].size();
165 rightcols = (msize_t)right[0].size();
166 if (leftcols && rightcols && (rightrows == leftcols))
167 {
168 resize(leftrows, rightcols, (*tresult), true);//could use zeros(), but common index last lets us zero at the same time
169 msize_t i, j, k;
170 for (i = 0; i < leftrows; ++i)
171 {
172 for (j = 0; j < rightcols; ++j)
173 {
174 A accum = 0;
175 for (k = 0; k < leftcols; ++k)
176 {
177 accum += left[i][k] * right[k][j];
178 }
179 (*tresult)[i][j] = accum;
180 }
181 }
182 } else {
183 result.resize(0);
184 return;
185 }
186 } else {
187 result.resize(0);
188 return;
189 }
190 if (copyout)
191 {
192 result = tempstorage;
193 }
194 }
195
196 template <typename T1, typename T2, typename T3>
197 void MatrixFunctions::multiply(const vector<vector<T1> > &left, const T2 right, vector<vector<T3> > &result)
198 {
199 msize_t leftrows = (msize_t)left.size(), leftcols;
200 bool doresize = true;
201 if (&left == &result)
202 {
203 doresize = false;//don't resize if an input is an output
204 }
205 if (leftrows)
206 {
207 leftcols = (msize_t)left[0].size();
208 if (leftcols)
209 {
210 if (doresize) resize(leftrows, leftcols, result, true);
211 msize_t i, j;
212 for (i = 0; i < leftrows; ++i)
213 {
214 for (j = 0; j < leftcols; ++j)
215 {
216 result[i][j] = left[i][j] * right;
217 }
218 }
219 } else {
220 result.resize(0);
221 return;
222 }
223 } else {
224 result.resize(0);
225 return;
226 }
227 }
228
229 template<typename T>
230 void MatrixFunctions::rref_big(vector<vector<T> > &inout)
231 {
232 msize_t rows = (msize_t)inout.size(), cols;
233 if (rows > 0)
234 {
235 cols = (msize_t)inout[0].size();
236 if (cols > 0)
237 {
238 vector<msize_t> pivots(rows, -1), missingPivots;
239 msize_t i, j, k, myrow = 0;
240 msize_t pivotrow;
241 T tempval;
242 for (i = 0; i < cols; ++i)
243 {
244 if (myrow >= rows) break;//no pivots left
245 tempval = 0;
246 pivotrow = -1;
247 for (j = myrow; j < rows; ++j)
248 {//only search below for new pivot
249 if (abs(inout[j][i]) > tempval)
250 {
251 pivotrow = (msize_t)j;
252 tempval = abs(inout[j][i]);
253 }
254 }
255 if (pivotrow == -1)
256 {//naively expect linearly dependence to show as an exact zero
257 missingPivots.push_back(i);//record the missing pivot
258 continue;//move to the next column
259 }
260 inout[pivotrow].swap(inout[myrow]);//STL swap via pointers for constant time row swap
261 pivots[myrow] = i;//save the pivot location for back substitution
262 tempval = inout[myrow][i];
263 inout[myrow][i] = (T)1;
264 for (j = i + 1; j < cols; ++j)
265 {
266 inout[myrow][j] /= tempval;//divide row by pivot
267 }
268 for (j = myrow + 1; j < rows; ++j)
269 {//zero ONLY below pivot for now
270 tempval = inout[j][i];
271 inout[j][i] = (T)0;
272 for (k = i + 1; k < cols; ++k)
273 {
274 inout[j][k] -= tempval * inout[myrow][k];
275 }
276 }
277 ++myrow;//increment row on successful pivot
278 }
279 msize_t numMissing = (msize_t)missingPivots.size();
280 if (myrow > 1)//if there is only 1 pivot, there is no back substitution to do
281 {
282 msize_t lastPivotCol = pivots[myrow - 1];
283 for (i = myrow - 1; i > 0; --i)//loop through pivots, can't zero above the top pivot so exclude it
284 {
285 msize_t pivotCol = pivots[i];
286 for (j = i - 1; j >= 0; --j)//loop through rows above pivot
287 {
288 tempval = inout[j][pivotCol];
289 inout[j][pivotCol] = (T)0;//flat zero the entry above the pivot
290 for (k = numMissing - 1; k >= 0; --k)//back substitute within pivot range where pivots are missing
291 {
292 msize_t missingCol = missingPivots[k];
293 if (missingCol <= pivotCol) break;//equals will never trip, but whatever
294 inout[j][missingCol] -= tempval * inout[i][missingCol];
295 }
296 for (k = lastPivotCol + 1; k < cols; ++k)//loop through elements that are outside the pivot area
297 {
298 inout[j][k] -= tempval * inout[i][k];
299 }
300 }
301 }
302 }
303 } else {
304 inout.resize(0);
305 return;
306 }
307 } else {
308 inout.resize(0);
309 return;
310 }
311 }
312
313 template<typename T>
314 void MatrixFunctions::rref(vector<vector<T> > &inout)
315 {
316 msize_t rows = (msize_t)inout.size(), cols;
317 if (rows)
318 {
319 cols = (msize_t)inout[0].size();
320 if (cols)
321 {
322 if (rows > 7 || cols > 7)//when the matrix has this many rows/columns, it is faster to allocate storage for tracking pivots, and back substitute
323 {
324 rref_big(inout);
325 return;
326 }
327 msize_t i, j, k, myrow = 0;
328 msize_t pivotrow;
329 T tempval;
330 for (i = 0; i < cols; ++i)
331 {
332 if (myrow >= rows) break;//no pivots left
333 tempval = 0;
334 pivotrow = -1;
335 for (j = myrow; j < rows; ++j)
336 {//only search below for new pivot
337 if (abs(inout[j][i]) > tempval)
338 {
339 pivotrow = (msize_t)j;
340 tempval = abs(inout[j][i]);
341 }
342 }
343 if (pivotrow == -1)//it may be a good idea to include a "very small value" check here, but it could mess up if used on a matrix with all values very small
344 {//naively expect linearly dependence to show as an exact zero
345 continue;//move to the next column
346 }
347 inout[pivotrow].swap(inout[myrow]);//STL swap via pointers for constant time row swap
348 tempval = inout[myrow][i];
349 inout[myrow][i] = 1;
350 for (j = i + 1; j < cols; ++j)
351 {
352 inout[myrow][j] /= tempval;//divide row by pivot
353 }
354 for (j = 0; j < myrow; ++j)
355 {//zero above pivot
356 tempval = inout[j][i];
357 inout[j][i] = 0;
358 for (k = i + 1; k < cols; ++k)
359 {
360 inout[j][k] -= tempval * inout[myrow][k];
361 }
362 }
363 for (j = myrow + 1; j < rows; ++j)
364 {//zero below pivot
365 tempval = inout[j][i];
366 inout[j][i] = 0;
367 for (k = i + 1; k < cols; ++k)
368 {
369 inout[j][k] -= tempval * inout[myrow][k];
370 }
371 }
372 ++myrow;//increment row on successful pivot
373 }
374 } else {
375 inout.resize(0);
376 return;
377 }
378 } else {
379 inout.resize(0);
380 return;
381 }
382 }
383
384 template<typename T>
385 void MatrixFunctions::inverse(const vector<vector<T> > &in, vector<vector<T> > &result)
386 {//rref implementation, there are faster (more complicated) ways - if it isn't invertible, it will hand back something strange
387 msize_t inrows = (msize_t)in.size(), incols;
388 if (inrows)
389 {
390 incols = (msize_t)in[0].size();
391 if (incols == inrows)
392 {
393 vector<vector<T> > inter, inter2;
394 identity(incols, inter2);
395 horizCat(in, inter2, inter);
396 rref(inter);
397 getChunk(0, inrows, incols, incols * 2, inter, result);//already using a local variable, doesn't need to check for reference duplicity
398 } else {
399 result.resize(0);
400 return;
401 }
402 } else {
403 result.resize(0);
404 return;
405 }
406 }
407
408 template <typename T1, typename T2, typename T3>
409 void MatrixFunctions::add(const vector<vector<T1> >& left, const vector<vector<T2> >& right, vector<vector<T3> >& result)
410 {
411 msize_t inrows = (msize_t)left.size(), incols;
412 bool doresize = true;
413 if (&left == &result || &right == &result)
414 {
415 doresize = false;//don't resize if an input is an output - this is ok for addition, don't need a copy
416 }
417 if (inrows)
418 {
419 incols = (msize_t)left[0].size();
420 if (inrows == (msize_t)right.size() && incols == (msize_t)right[0].size())
421 {
422 if (doresize) resize(inrows, incols, result, true);
423 for (msize_t i = 0; i < inrows; ++i)
424 {
425 for (msize_t j = 0; j < incols; ++j)
426 {
427 result[i][j] = left[i][j] + right[i][j];
428 }
429 }
430 } else {
431 result.resize(0);//use empty matrix for error condition
432 return;
433 }
434 } else {
435 result.resize(0);
436 return;
437 }
438 }
439
440 template <typename T1, typename T2, typename T3>
441 void MatrixFunctions::add(const vector<vector<T1> >& left, const T2 right, vector<vector<T3> >& result)
442 {
443 msize_t inrows = (msize_t)left.size(), incols;
444 bool doresize = true;
445 if (&left == &result)
446 {
447 doresize = false;//don't resize if an input is an output - this is ok for addition, don't need a copy
448 }
449 if (inrows)
450 {
451 incols = (msize_t)left[0].size();
452 if (doresize) resize(inrows, incols, result, true);
453 for (msize_t i = 0; i < inrows; ++i)
454 {
455 for (msize_t j = 0; j < incols; ++j)
456 {
457 result[i][j] = left[i][j] + right;
458 }
459 }
460 } else {
461 result.resize(0);
462 return;
463 }
464 }
465
466 template <typename T1, typename T2, typename T3>
467 void MatrixFunctions::subtract(const vector<vector<T1> >& left, const vector<vector<T2> >& right, vector<vector<T3> >& result)
468 {
469 msize_t inrows = (msize_t)left.size(), incols;
470 bool doresize = true;
471 if (&left == &result || &right == &result)
472 {
473 doresize = false;//don't resize if an input is an output
474 }
475 if (inrows)
476 {
477 incols = (msize_t)left[0].size();
478 if (inrows == (msize_t)right.size() && incols == (msize_t)right[0].size())
479 {
480 if (doresize) resize(inrows, incols, result, true);
481 for (msize_t i = 0; i < inrows; ++i)
482 {
483 for (msize_t j = 0; j < incols; ++j)
484 {
485 result[i][j] = left[i][j] - right[i][j];
486 }
487 }
488 } else {
489 result.resize(0);
490 return;
491 }
492 } else {
493 result.resize(0);
494 return;
495 }
496 }
497
498 template<typename T>
499 void MatrixFunctions::transpose(const vector<vector<T> > &in, vector<vector<T> > &result)
500 {
501 msize_t inrows = (msize_t)in.size(), incols;
502 vector<vector<T> > tempstorage, *tresult = &result;
503 bool copyout = false;
504 if (&in == &result)
505 {
506 copyout = true;
507 tresult = &tempstorage;
508 }
509 if (inrows)
510 {
511 incols = (msize_t)in[0].size();
512 resize(incols, inrows, (*tresult), true);
513 for (msize_t i = 0; i < inrows; ++i)
514 {
515 for (msize_t j = 0; j < incols; ++j)
516 {
517 (*tresult)[j][i] = in[i][j];
518 }
519 }
520 } else {
521 result.resize(0);
522 }
523 if (copyout)
524 {
525 result = tempstorage;
526 }
527 }
528
529 template<typename T>
530 bool MatrixFunctions::checkDim(const vector<vector<T> > &in)
531 {
532 bool ret = true;
533 msize_t rows = (msize_t)in.size(), columns;
534 if (rows)
535 {
536 columns = (msize_t)in[0].size();
537 for (msize_t i = 1; i < rows; ++i)
538 {
539 if (in[i].size() != columns)
540 {
541 ret = false;
542 }
543 }
544 }
545 return ret;
546 }
547
548 template<typename T>
549 void MatrixFunctions::resize(const msize_t rows, const msize_t columns, vector<vector<T> >& result, bool destructive)
550 {
551 if (destructive && result.size() && ((msize_t)result.capacity() < rows || (msize_t)result[0].capacity() < columns))
552 {//for large matrices, copying to preserve contents is slow
553 result.resize(0);//not intended to dealloc, just to set number of items to copy to zero
554 }//default is nondestructive resize, copies everything
555 result.resize(rows);
556 for (msize_t i = 0; i < (const msize_t)rows; ++i)
557 {//naive method, may end up copying everything twice if both row and col resizes require realloc
558 result[i].resize(columns);
559 }
560 }
561
562 template<typename T>
563 void MatrixFunctions::zeros(const msize_t rows, const msize_t columns, vector<vector<T> >& result)
564 {
565 resize(rows, columns, result, true);
566 for (msize_t i = 0; i < rows; ++i)
567 {
568 for (msize_t j = 0; j < columns; ++j)
569 {
570 result[i][j] = 0;//should cast to float or double fine
571 }
572 }
573 }
574
575 template<typename T>
576 void MatrixFunctions::ones(const msize_t rows, const msize_t columns, vector<vector<T> >& result)
577 {
578 resize(rows, columns, result, true);
579 for (msize_t i = 0; i < rows; ++i)
580 {
581 for (msize_t j = 0; j < columns; ++j)
582 {
583 result[i][j] = 1;//should cast to float or double fine
584 }
585 }
586 }
587
588 template<typename T>
589 void MatrixFunctions::identity(const msize_t size, vector<vector<T> >& result)
590 {
591 resize(size, size, result, true);
592 for (msize_t i = 0; i < (const msize_t)size; ++i)
593 {
594 for (msize_t j = 0; j < (const msize_t)size; ++j)
595 {
596 result[i][j] = ((i == j) ? 1 : 0);//ditto, forgive the ternary
597 }
598 }
599 }
600
601 template <typename T1, typename T2, typename T3>
602 void MatrixFunctions::horizCat(const vector<vector<T1> >& left, const vector<vector<T2> >& right, vector<vector<T3> >& result)
603 {
604 msize_t inrows = (msize_t)left.size(), leftcols, rightcols;
605 vector<vector<T3> > tempstorage, *tresult = &result;
606 bool copyout = false;
607 if (&left == &result || &right == &result)
608 {
609 copyout = true;
610 tresult = &tempstorage;
611 }
612 if (inrows && inrows == (msize_t)right.size())
613 {
614 leftcols = (msize_t)left[0].size();
615 rightcols = (msize_t)right[0].size();
616 (*tresult) = left;//use STL copy to start
617 resize(inrows, leftcols + rightcols, (*tresult));//values survive nondestructive resize
618 for (msize_t i = 0; i < inrows; ++i)
619 {
620 for (msize_t j = 0; j < rightcols; ++j)
621 {
622 (*tresult)[i][j + leftcols] = right[i][j];
623 }
624 }
625 } else {
626 result.resize(0);
627 return;
628 }
629 if (copyout)
630 {
631 result = tempstorage;
632 }
633 }
634
635 template <typename T1, typename T2, typename T3>
636 void MatrixFunctions::vertCat(const vector<vector<T1> >& top, const vector<vector<T2> >& bottom, vector<vector<T3> >& result)
637 {
638 msize_t toprows = (msize_t)top.size(), botrows = (msize_t)bottom.size(), incols;
639 vector<vector<T3> > tempstorage, *tresult = &result;
640 bool copyout = false;
641 if (&top == &result || &bottom == &result)
642 {
643 copyout = true;
644 tresult = &tempstorage;
645 }
646 if (toprows && botrows)
647 {
648 incols = (msize_t)top[0].size();
649 if (incols == (msize_t)bottom[0].size())
650 {
651 (*tresult) = top;
652 resize(toprows + botrows, incols, (*tresult));//nondestructive resize
653 for (msize_t i = 0; i < botrows; ++i)
654 {
655 for (msize_t j = 0; j < incols; ++j)
656 {
657 (*tresult)[i + toprows][j] = bottom[i][j];
658 }
659 }
660 } else {
661 result.resize(0);
662 return;
663 }
664 } else {
665 result.resize(0);
666 return;
667 }
668 if (copyout)
669 {
670 result = tempstorage;
671 }
672 }
673
674 template<typename T>
675 void MatrixFunctions::getChunk(const msize_t firstrow, const msize_t lastrow, const msize_t firstcol, const msize_t lastcol, const vector<vector<T> >& in, vector<vector<T> >& result)
676 {
677 msize_t outrows = lastrow - firstrow;
678 msize_t outcols = lastcol - firstcol;
679 if (lastrow <= firstrow || lastcol <= firstcol || firstrow < 0 || firstcol < 0 || lastrow > (msize_t)in.size() || lastcol > (msize_t)in[0].size())
680 {
681 result.resize(0);
682 return;
683 }
684 vector<vector<T> > tempstorage, *tresult = &result;
685 bool copyout = false;
686 if (&in == &result)
687 {
688 copyout = true;
689 tresult = &tempstorage;
690 }
691 resize(outrows, outcols, (*tresult), true);
692 for (msize_t i = 0; i < outrows; ++i)
693 {
694 for (msize_t j = 0; j < outcols; ++j)
695 {
696 (*tresult)[i][j] = in[i + firstrow][j + firstcol];
697 }
698 }
699 if (copyout)
700 {
701 result = tempstorage;
702 }
703 }
704
705}
706
707#endif
708
Definition MatrixFunctions.h:43
static void resize(const msize_t rows, const msize_t columns, vector< vector< T > > &result, bool destructive=false)
Definition MatrixFunctions.h:549
static void vertCat(const vector< vector< T1 > > &top, const vector< vector< T2 > > &bottom, vector< vector< T3 > > &result)
Definition MatrixFunctions.h:636
static void rref(vector< vector< T > > &inout)
Definition MatrixFunctions.h:314
static void ones(const msize_t rows, const msize_t columns, vector< vector< T > > &result)
Definition MatrixFunctions.h:576
static void add(const vector< vector< T1 > > &left, const vector< vector< T2 > > &right, vector< vector< T3 > > &result)
Definition MatrixFunctions.h:409
static void transpose(const vector< vector< T > > &in, vector< vector< T > > &result)
Definition MatrixFunctions.h:499
static void zeros(const msize_t rows, const msize_t columns, vector< vector< T > > &result)
Definition MatrixFunctions.h:563
static void identity(const msize_t size, vector< vector< T > > &result)
Definition MatrixFunctions.h:589
static void multiply(const vector< vector< T1 > > &left, const vector< vector< T2 > > &right, vector< vector< T3 > > &result)
Definition MatrixFunctions.h:152
static bool checkDim(const vector< vector< T > > &in)
Definition MatrixFunctions.h:530
static void subtract(const vector< vector< T1 > > &left, const vector< vector< T2 > > &right, vector< vector< T3 > > &result)
Definition MatrixFunctions.h:467
static void horizCat(const vector< vector< T1 > > &left, const vector< vector< T2 > > &right, vector< vector< T3 > > &result)
Definition MatrixFunctions.h:602
static void getChunk(const msize_t firstrow, const msize_t lastrow, const msize_t firstcol, const msize_t lastcol, const vector< vector< T > > &in, vector< vector< T > > &result)
Definition MatrixFunctions.h:675
static void inverse(const vector< vector< T > > &in, vector< vector< T > > &result)
Definition MatrixFunctions.h:385
namespace for all CiftiLib functionality
Definition CiftiBrainModelsMap.h:42