ARTS 2.5.4 (git: 4c0d3b4d)
matpackIV.cc
Go to the documentation of this file.
1/* Copyright (C) 2002-2012 Stefan Buehler <sbuehler@ltu.se>
2
3 This program is free software; you can redistribute it and/or modify it
4 under the terms of the GNU General Public License as published by the
5 Free Software Foundation; either version 2, or (at your option) any
6 later version.
7
8 This program is distributed in the hope that it will be useful,
9 but WITHOUT ANY WARRANTY; without even the implied warranty of
10 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11 GNU General Public License for more details.
12
13 You should have received a copy of the GNU General Public License
14 along with this program; if not, write to the Free Software
15 Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
16 USA. */
17
25#include "matpackIV.h"
26
27#include "exceptions.h"
28
31Tensor3View* Iterator4D::operator->() { return &msv; }
32
34Tensor3View& Iterator4D::operator*() { return msv; }
35
38const ConstTensor3View* ConstIterator4D::operator->() const { return &msv; }
39
41const ConstTensor3View& ConstIterator4D::operator*() const { return msv; }
42
43// Functions for ConstTensor4View:
44// ------------------------------
45
49ConstTensor4View ConstTensor4View::operator()(const Range& b,
50 const Range& p,
51 const Range& r,
52 const Range& c) const {
53 return ConstTensor4View(mdata, mbr, mpr, mrr, mcr, b, p, r, c);
54}
55
58ConstTensor3View ConstTensor4View::operator()(const Range& b,
59 const Range& p,
60 const Range& r,
61 Index c) const {
62 // Check that c is valid:
63 ARTS_ASSERT(0 <= c);
64 ARTS_ASSERT(c < mcr.mextent);
65
66 return ConstTensor3View(
67 mdata + mcr.mstart + c * mcr.mstride, mbr, mpr, mrr, b, p, r);
68}
69
72ConstTensor3View ConstTensor4View::operator()(const Range& b,
73 const Range& p,
74 Index r,
75 const Range& c) const {
76 // Check that r is valid:
77 ARTS_ASSERT(0 <= r);
78 ARTS_ASSERT(r < mrr.mextent);
79
80 return ConstTensor3View(
81 mdata + mrr.mstart + r * mrr.mstride, mbr, mpr, mcr, b, p, c);
82}
83
86ConstTensor3View ConstTensor4View::operator()(const Range& b,
87 Index p,
88 const Range& r,
89 const Range& c) const {
90 // Check that p is valid:
91 ARTS_ASSERT(0 <= p);
92 ARTS_ASSERT(p < mpr.mextent);
93
94 return ConstTensor3View(
95 mdata + mpr.mstart + p * mpr.mstride, mbr, mrr, mcr, b, r, c);
96}
97
100ConstTensor3View ConstTensor4View::operator()(Index b,
101 const Range& p,
102 const Range& r,
103 const Range& c) const {
104 // Check that b is valid:
105 ARTS_ASSERT(0 <= b);
106 ARTS_ASSERT(b < mbr.mextent);
107
108 return ConstTensor3View(
109 mdata + mbr.mstart + b * mbr.mstride, mpr, mrr, mcr, p, r, c);
110}
111
114ConstMatrixView ConstTensor4View::operator()(const Range& b,
115 const Range& p,
116 Index r,
117 Index c) const {
118 // Check that r and c are valid:
119 ARTS_ASSERT(0 <= r);
120 ARTS_ASSERT(0 <= c);
121 ARTS_ASSERT(r < mrr.mextent);
122 ARTS_ASSERT(c < mcr.mextent);
123
124 return ConstMatrixView(
125 mdata + mrr.mstart + r * mrr.mstride + mcr.mstart + c * mcr.mstride,
126 mbr,
127 mpr,
128 b,
129 p);
130}
131
134ConstMatrixView ConstTensor4View::operator()(const Range& b,
135 Index p,
136 const Range& r,
137 Index c) const {
138 // Check that p and c are valid:
139 ARTS_ASSERT(0 <= p);
140 ARTS_ASSERT(0 <= c);
141 ARTS_ASSERT(p < mpr.mextent);
142 ARTS_ASSERT(c < mcr.mextent);
143
144 return ConstMatrixView(
145 mdata + mpr.mstart + p * mpr.mstride + mcr.mstart + c * mcr.mstride,
146 mbr,
147 mrr,
148 b,
149 r);
150}
151
154ConstMatrixView ConstTensor4View::operator()(const Range& b,
155 Index p,
156 Index r,
157 const Range& c) const {
158 // Check that p and r are valid:
159 ARTS_ASSERT(0 <= p);
160 ARTS_ASSERT(0 <= r);
161 ARTS_ASSERT(p < mpr.mextent);
162 ARTS_ASSERT(r < mrr.mextent);
163
164 return ConstMatrixView(
165 mdata + mpr.mstart + p * mpr.mstride + mrr.mstart + r * mrr.mstride,
166 mbr,
167 mcr,
168 b,
169 c);
170}
171
174ConstMatrixView ConstTensor4View::operator()(Index b,
175 const Range& p,
176 Index r,
177 const Range& c) const {
178 // Check that b and r are valid:
179 ARTS_ASSERT(0 <= b);
180 ARTS_ASSERT(0 <= r);
181 ARTS_ASSERT(b < mbr.mextent);
182 ARTS_ASSERT(r < mrr.mextent);
183
184 return ConstMatrixView(
185 mdata + mbr.mstart + b * mbr.mstride + mrr.mstart + r * mrr.mstride,
186 mpr,
187 mcr,
188 p,
189 c);
190}
191
194ConstMatrixView ConstTensor4View::operator()(Index b,
195 const Range& p,
196 const Range& r,
197 Index c) const {
198 // Check that b and c are valid:
199 ARTS_ASSERT(0 <= b);
200 ARTS_ASSERT(0 <= c);
201 ARTS_ASSERT(b < mbr.mextent);
202 ARTS_ASSERT(c < mcr.mextent);
203
204 return ConstMatrixView(
205 mdata + mbr.mstart + b * mbr.mstride + mcr.mstart + c * mcr.mstride,
206 mpr,
207 mrr,
208 p,
209 r);
210}
211
214ConstMatrixView ConstTensor4View::operator()(Index b,
215 Index p,
216 const Range& r,
217 const Range& c) const {
218 // Check that b and p are valid:
219 ARTS_ASSERT(0 <= b);
220 ARTS_ASSERT(0 <= p);
221 ARTS_ASSERT(b < mbr.mextent);
222 ARTS_ASSERT(p < mpr.mextent);
223
224 return ConstMatrixView(
225 mdata + mbr.mstart + b * mbr.mstride + mpr.mstart + p * mpr.mstride,
226 mrr,
227 mcr,
228 r,
229 c);
230}
231
234ConstVectorView ConstTensor4View::operator()(const Range& b,
235 Index p,
236 Index r,
237 Index c) const {
238 // Check that p, r and c are valid:
239 ARTS_ASSERT(0 <= p);
240 ARTS_ASSERT(0 <= r);
241 ARTS_ASSERT(0 <= c);
242 ARTS_ASSERT(p < mpr.mextent);
243 ARTS_ASSERT(r < mrr.mextent);
244 ARTS_ASSERT(c < mcr.mextent);
245
246 return ConstVectorView(mdata + mpr.mstart + p * mpr.mstride + mrr.mstart +
247 r * mrr.mstride + mcr.mstart + c * mcr.mstride,
248 mbr,
249 b);
250}
251
254ConstVectorView ConstTensor4View::operator()(Index b,
255 const Range& p,
256 Index r,
257 Index c) const {
258 // Check that b, r and c are valid:
259 ARTS_ASSERT(0 <= b);
260 ARTS_ASSERT(0 <= r);
261 ARTS_ASSERT(0 <= c);
262 ARTS_ASSERT(b < mbr.mextent);
263 ARTS_ASSERT(r < mrr.mextent);
264 ARTS_ASSERT(c < mcr.mextent);
265
266 return ConstVectorView(mdata + mbr.mstart + b * mbr.mstride + mrr.mstart +
267 r * mrr.mstride + mcr.mstart + c * mcr.mstride,
268 mpr,
269 p);
270}
271
274ConstVectorView ConstTensor4View::operator()(Index b,
275 Index p,
276 const Range& r,
277 Index c) const {
278 // Check that b, p and c are valid:
279 ARTS_ASSERT(0 <= b);
280 ARTS_ASSERT(0 <= p);
281 ARTS_ASSERT(0 <= c);
282 ARTS_ASSERT(b < mbr.mextent);
283 ARTS_ASSERT(p < mpr.mextent);
284 ARTS_ASSERT(c < mcr.mextent);
285
286 return ConstVectorView(mdata + mbr.mstart + b * mbr.mstride + mpr.mstart +
287 p * mpr.mstride + mcr.mstart + c * mcr.mstride,
288 mrr,
289 r);
290}
291
294ConstVectorView ConstTensor4View::operator()(Index b,
295 Index p,
296 Index r,
297 const Range& c) const {
298 // Check that b, p and r are valid:
299 ARTS_ASSERT(0 <= b);
300 ARTS_ASSERT(0 <= p);
301 ARTS_ASSERT(0 <= r);
302 ARTS_ASSERT(b < mbr.mextent);
303 ARTS_ASSERT(p < mpr.mextent);
304 ARTS_ASSERT(r < mrr.mextent);
305
306 return ConstVectorView(mdata + mbr.mstart + b * mbr.mstride + mpr.mstart +
307 p * mpr.mstride + mrr.mstart + r * mrr.mstride,
308 mcr,
309 c);
310}
311
318Numeric* Tensor4View::get_c_array() ARTS_NOEXCEPT {
319 ARTS_ASSERT(mbr.mstart == 0 and
320 (mbr.mstride == mcr.mextent * mrr.mextent * mpr.mextent or
321 size() == 0),
322 "Book ",
323 mbr)
324 ARTS_ASSERT(mpr.mstart == 0 and
325 (mpr.mstride == mcr.mextent * mrr.mextent or size() == 0),
326 "Page ",
327 mpr)
328 ARTS_ASSERT(mrr.mstart == 0 and (mrr.mstride == mcr.mextent or size() == 0),
329 "Row ",
330 mrr)
331 ARTS_ASSERT(
332 mcr.mstart == 0 and (mcr.mstride == 1 or size() == 0), "Column ", mcr)
333 return mdata;
334}
335
342const Numeric* Tensor4View::get_c_array() const ARTS_NOEXCEPT {
343 ARTS_ASSERT(mbr.mstart == 0 and
344 (mbr.mstride == mcr.mextent * mrr.mextent * mpr.mextent or
345 size() == 0),
346 "Book ",
347 mbr)
348 ARTS_ASSERT(mpr.mstart == 0 and
349 (mpr.mstride == mcr.mextent * mrr.mextent or size() == 0),
350 "Page ",
351 mpr)
352 ARTS_ASSERT(mrr.mstart == 0 and (mrr.mstride == mcr.mextent or size() == 0),
353 "Row ",
354 mrr)
355 ARTS_ASSERT(
356 mcr.mstart == 0 and (mcr.mstride == 1 or size() == 0), "Column ", mcr)
357 return mdata;
358}
359
361ConstIterator4D ConstTensor4View::begin() const {
362 return ConstIterator4D(ConstTensor3View(mdata + mbr.mstart, mpr, mrr, mcr),
363 mbr.mstride);
364}
365
367ConstIterator4D ConstTensor4View::end() const {
368 return ConstIterator4D(
369 ConstTensor3View(
370 mdata + mbr.mstart + (mbr.mextent) * mbr.mstride, mpr, mrr, mcr),
371 mbr.mstride);
372}
373
375ConstTensor4View::ConstTensor4View(const ConstTensor3View& a)
376 : mbr(0, 1, a.mpr.mextent * a.mrr.mextent * a.mcr.mextent),
377 mpr(a.mpr),
378 mrr(a.mrr),
379 mcr(a.mcr),
380 mdata(a.mdata) {
381 // Nothing to do here.
382}
383
388ConstTensor4View::ConstTensor4View(Numeric* data,
389 const Range& br,
390 const Range& pr,
391 const Range& rr,
392 const Range& cr)
393 : mbr(br), mpr(pr), mrr(rr), mcr(cr), mdata(data) {
394 // Nothing to do here.
395}
396
404ConstTensor4View::ConstTensor4View(Numeric* data,
405 const Range& pb,
406 const Range& pp,
407 const Range& pr,
408 const Range& pc,
409 const Range& nb,
410 const Range& np,
411 const Range& nr,
412 const Range& nc)
413 : mbr(pb, nb), mpr(pp, np), mrr(pr, nr), mcr(pc, nc), mdata(data) {
414 // Nothing to do here.
415}
416
420std::ostream& operator<<(std::ostream& os, const ConstTensor4View& v) {
421 // Page iterators:
422 ConstIterator4D ib = v.begin();
423 const ConstIterator4D end_book = v.end();
424
425 if (ib != end_book) {
426 os << *ib;
427 ++ib;
428 }
429
430 for (; ib != end_book; ++ib) {
431 os << "\n\n";
432 os << *ib;
433 }
434
435 return os;
436}
437
438// Functions for Tensor4View:
439// -------------------------
440
444Tensor4View Tensor4View::operator()(const Range& b,
445 const Range& p,
446 const Range& r,
447 const Range& c) {
448 return Tensor4View(mdata, mbr, mpr, mrr, mcr, b, p, r, c);
449}
450
453Tensor3View Tensor4View::operator()(const Range& b,
454 const Range& p,
455 const Range& r,
456 Index c) {
457 // Check that c is valid:
458 ARTS_ASSERT(0 <= c);
459 ARTS_ASSERT(c < mcr.mextent);
460
461 return Tensor3View(
462 mdata + mcr.mstart + c * mcr.mstride, mbr, mpr, mrr, b, p, r);
463}
464
467Tensor3View Tensor4View::operator()(const Range& b,
468 const Range& p,
469 Index r,
470 const Range& c) {
471 // Check that r is valid:
472 ARTS_ASSERT(0 <= r);
473 ARTS_ASSERT(r < mrr.mextent);
474
475 return Tensor3View(
476 mdata + mrr.mstart + r * mrr.mstride, mbr, mpr, mcr, b, p, c);
477}
478
481Tensor3View Tensor4View::operator()(const Range& b,
482 Index p,
483 const Range& r,
484 const Range& c) {
485 // Check that p is valid:
486 ARTS_ASSERT(0 <= p);
487 ARTS_ASSERT(p < mpr.mextent);
488
489 return Tensor3View(
490 mdata + mpr.mstart + p * mpr.mstride, mbr, mrr, mcr, b, r, c);
491}
492
495Tensor3View Tensor4View::operator()(Index b,
496 const Range& p,
497 const Range& r,
498 const Range& c) {
499 // Check that b is valid:
500 ARTS_ASSERT(0 <= b);
501 ARTS_ASSERT(b < mbr.mextent);
502
503 return Tensor3View(
504 mdata + mbr.mstart + b * mbr.mstride, mpr, mrr, mcr, p, r, c);
505}
506
509MatrixView Tensor4View::operator()(const Range& b,
510 const Range& p,
511 Index r,
512 Index c) {
513 // Check that r and c are valid:
514 ARTS_ASSERT(0 <= r);
515 ARTS_ASSERT(0 <= c);
516 ARTS_ASSERT(r < mrr.mextent);
517 ARTS_ASSERT(c < mcr.mextent);
518
519 return MatrixView(
520 mdata + mrr.mstart + r * mrr.mstride + mcr.mstart + c * mcr.mstride,
521 mbr,
522 mpr,
523 b,
524 p);
525}
526
529MatrixView Tensor4View::operator()(const Range& b,
530 Index p,
531 const Range& r,
532 Index c) {
533 // Check that p and c are valid:
534 ARTS_ASSERT(0 <= p);
535 ARTS_ASSERT(0 <= c);
536 ARTS_ASSERT(p < mpr.mextent);
537 ARTS_ASSERT(c < mcr.mextent);
538
539 return MatrixView(
540 mdata + mpr.mstart + p * mpr.mstride + mcr.mstart + c * mcr.mstride,
541 mbr,
542 mrr,
543 b,
544 r);
545}
546
549MatrixView Tensor4View::operator()(const Range& b,
550 Index p,
551 Index r,
552 const Range& c) {
553 // Check that p and r are valid:
554 ARTS_ASSERT(0 <= p);
555 ARTS_ASSERT(0 <= r);
556 ARTS_ASSERT(p < mpr.mextent);
557 ARTS_ASSERT(r < mrr.mextent);
558
559 return MatrixView(
560 mdata + mpr.mstart + p * mpr.mstride + mrr.mstart + r * mrr.mstride,
561 mbr,
562 mcr,
563 b,
564 c);
565}
566
569MatrixView Tensor4View::operator()(Index b,
570 const Range& p,
571 Index r,
572 const Range& c) {
573 // Check that b and r are valid:
574 ARTS_ASSERT(0 <= b);
575 ARTS_ASSERT(0 <= r);
576 ARTS_ASSERT(b < mbr.mextent);
577 ARTS_ASSERT(r < mrr.mextent);
578
579 return MatrixView(
580 mdata + mbr.mstart + b * mbr.mstride + mrr.mstart + r * mrr.mstride,
581 mpr,
582 mcr,
583 p,
584 c);
585}
586
589MatrixView Tensor4View::operator()(Index b,
590 const Range& p,
591 const Range& r,
592 Index c) {
593 // Check that b and c are valid:
594 ARTS_ASSERT(0 <= b);
595 ARTS_ASSERT(0 <= c);
596 ARTS_ASSERT(b < mbr.mextent);
597 ARTS_ASSERT(c < mcr.mextent);
598
599 return MatrixView(
600 mdata + mbr.mstart + b * mbr.mstride + mcr.mstart + c * mcr.mstride,
601 mpr,
602 mrr,
603 p,
604 r);
605}
606
609MatrixView Tensor4View::operator()(Index b,
610 Index p,
611 const Range& r,
612 const Range& c) {
613 // Check that b and p are valid:
614 ARTS_ASSERT(0 <= b);
615 ARTS_ASSERT(0 <= p);
616 ARTS_ASSERT(b < mbr.mextent);
617 ARTS_ASSERT(p < mpr.mextent);
618
619 return MatrixView(
620 mdata + mbr.mstart + b * mbr.mstride + mpr.mstart + p * mpr.mstride,
621 mrr,
622 mcr,
623 r,
624 c);
625}
626
629VectorView Tensor4View::operator()(const Range& b, Index p, Index r, Index c) {
630 // Check that p, r and c are valid:
631 ARTS_ASSERT(0 <= p);
632 ARTS_ASSERT(0 <= r);
633 ARTS_ASSERT(0 <= c);
634 ARTS_ASSERT(p < mpr.mextent);
635 ARTS_ASSERT(r < mrr.mextent);
636 ARTS_ASSERT(c < mcr.mextent);
637
638 return VectorView(mdata + mpr.mstart + p * mpr.mstride + mrr.mstart +
639 r * mrr.mstride + mcr.mstart + c * mcr.mstride,
640 mbr,
641 b);
642}
643
646VectorView Tensor4View::operator()(Index b, const Range& p, Index r, Index c) {
647 // Check that b, r and c are valid:
648 ARTS_ASSERT(0 <= b);
649 ARTS_ASSERT(0 <= r);
650 ARTS_ASSERT(0 <= c);
651 ARTS_ASSERT(b < mbr.mextent);
652 ARTS_ASSERT(r < mrr.mextent);
653 ARTS_ASSERT(c < mcr.mextent);
654
655 return VectorView(mdata + mbr.mstart + b * mbr.mstride + mrr.mstart +
656 r * mrr.mstride + mcr.mstart + c * mcr.mstride,
657 mpr,
658 p);
659}
660
663VectorView Tensor4View::operator()(Index b, Index p, const Range& r, Index c) {
664 // Check that b, p and c are valid:
665 ARTS_ASSERT(0 <= b);
666 ARTS_ASSERT(0 <= p);
667 ARTS_ASSERT(0 <= c);
668 ARTS_ASSERT(b < mbr.mextent);
669 ARTS_ASSERT(p < mpr.mextent);
670 ARTS_ASSERT(c < mcr.mextent);
671
672 return VectorView(mdata + mbr.mstart + b * mbr.mstride + mpr.mstart +
673 p * mpr.mstride + mcr.mstart + c * mcr.mstride,
674 mrr,
675 r);
676}
677
680VectorView Tensor4View::operator()(Index b, Index p, Index r, const Range& c) {
681 // Check that b, p and r are valid:
682 ARTS_ASSERT(0 <= b);
683 ARTS_ASSERT(0 <= p);
684 ARTS_ASSERT(0 <= r);
685 ARTS_ASSERT(b < mbr.mextent);
686 ARTS_ASSERT(p < mpr.mextent);
687 ARTS_ASSERT(r < mrr.mextent);
688
689 return VectorView(mdata + mbr.mstart + b * mbr.mstride + mpr.mstart +
690 p * mpr.mstride + mrr.mstart + r * mrr.mstride,
691 mcr,
692 c);
693}
694
696Iterator4D Tensor4View::begin() {
697 return Iterator4D(Tensor3View(mdata + mbr.mstart, mpr, mrr, mcr),
698 mbr.mstride);
699}
700
702Iterator4D Tensor4View::end() {
703 return Iterator4D(
704 Tensor3View(
705 mdata + mbr.mstart + (mbr.mextent) * mbr.mstride, mpr, mrr, mcr),
706 mbr.mstride);
707}
708
713Tensor4View& Tensor4View::operator=(const ConstTensor4View& m) {
714 // Check that sizes are compatible:
715 ARTS_ASSERT(mbr.mextent == m.mbr.mextent);
716 ARTS_ASSERT(mpr.mextent == m.mpr.mextent);
717 ARTS_ASSERT(mrr.mextent == m.mrr.mextent);
718 ARTS_ASSERT(mcr.mextent == m.mcr.mextent);
719
720 copy(m.begin(), m.end(), begin());
721 return *this;
722}
723
729Tensor4View& Tensor4View::operator=(const Tensor4View& m) {
730 // Check that sizes are compatible:
731 ARTS_ASSERT(mbr.mextent == m.mbr.mextent);
732 ARTS_ASSERT(mpr.mextent == m.mpr.mextent);
733 ARTS_ASSERT(mrr.mextent == m.mrr.mextent);
734 ARTS_ASSERT(mcr.mextent == m.mcr.mextent);
735
736 copy(m.begin(), m.end(), begin());
737 return *this;
738}
739
743Tensor4View& Tensor4View::operator=(const Tensor4& m) {
744 // Check that sizes are compatible:
745 ARTS_ASSERT(mbr.mextent == m.mbr.mextent);
746 ARTS_ASSERT(mpr.mextent == m.mpr.mextent);
747 ARTS_ASSERT(mrr.mextent == m.mrr.mextent);
748 ARTS_ASSERT(mcr.mextent == m.mcr.mextent);
749
750 copy(m.begin(), m.end(), begin());
751 return *this;
752}
753
756Tensor4View& Tensor4View::operator=(Numeric x) {
757 copy(x, begin(), end());
758 return *this;
759}
760
761// Some little helper functions:
762//------------------------------
763
765Tensor4View& Tensor4View::operator*=(Numeric x) {
766 const Iterator4D eb = end();
767 for (Iterator4D b = begin(); b != eb; ++b) {
768 *b *= x;
769 }
770 return *this;
771}
772
774Tensor4View& Tensor4View::operator/=(Numeric x) {
775 const Iterator4D eb = end();
776 for (Iterator4D b = begin(); b != eb; ++b) {
777 *b /= x;
778 }
779 return *this;
780}
781
783Tensor4View& Tensor4View::operator+=(Numeric x) {
784 const Iterator4D eb = end();
785 for (Iterator4D b = begin(); b != eb; ++b) {
786 *b += x;
787 }
788 return *this;
789}
790
792Tensor4View& Tensor4View::operator-=(Numeric x) {
793 const Iterator4D eb = end();
794 for (Iterator4D b = begin(); b != eb; ++b) {
795 *b -= x;
796 }
797 return *this;
798}
799
801Tensor4View& Tensor4View::operator*=(const ConstTensor4View& x) {
802 ARTS_ASSERT(nbooks() == x.nbooks());
803 ARTS_ASSERT(npages() == x.npages());
804 ARTS_ASSERT(nrows() == x.nrows());
805 ARTS_ASSERT(ncols() == x.ncols());
806 ConstIterator4D xb = x.begin();
807 Iterator4D b = begin();
808 const Iterator4D eb = end();
809 for (; b != eb; ++b, ++xb) {
810 *b *= *xb;
811 }
812 return *this;
813}
814
816Tensor4View& Tensor4View::operator/=(const ConstTensor4View& x) {
817 ARTS_ASSERT(nbooks() == x.nbooks());
818 ARTS_ASSERT(npages() == x.npages());
819 ARTS_ASSERT(nrows() == x.nrows());
820 ARTS_ASSERT(ncols() == x.ncols());
821 ConstIterator4D xb = x.begin();
822 Iterator4D b = begin();
823 const Iterator4D eb = end();
824 for (; b != eb; ++b, ++xb) {
825 *b /= *xb;
826 }
827 return *this;
828}
829
831Tensor4View& Tensor4View::operator+=(const ConstTensor4View& x) {
832 ARTS_ASSERT(nbooks() == x.nbooks());
833 ARTS_ASSERT(npages() == x.npages());
834 ARTS_ASSERT(nrows() == x.nrows());
835 ARTS_ASSERT(ncols() == x.ncols());
836 ConstIterator4D xb = x.begin();
837 Iterator4D b = begin();
838 const Iterator4D eb = end();
839 for (; b != eb; ++b, ++xb) {
840 *b += *xb;
841 }
842 return *this;
843}
844
846Tensor4View& Tensor4View::operator-=(const ConstTensor4View& x) {
847 ARTS_ASSERT(nbooks() == x.nbooks());
848 ARTS_ASSERT(npages() == x.npages());
849 ARTS_ASSERT(nrows() == x.nrows());
850 ARTS_ASSERT(ncols() == x.ncols());
851 ConstIterator4D xb = x.begin();
852 Iterator4D b = begin();
853 const Iterator4D eb = end();
854 for (; b != eb; ++b, ++xb) {
855 *b -= *xb;
856 }
857 return *this;
858}
859
861Tensor4View::Tensor4View(const Tensor3View& a)
862 : ConstTensor4View(
863 a.mdata,
864 Range(0, 1, a.mpr.mextent * a.mrr.mextent * a.mcr.mextent),
865 a.mpr,
866 a.mrr,
867 a.mcr) {
868 // Nothing to do here.
869}
870
874Tensor4View::Tensor4View(Numeric* data,
875 const Range& br,
876 const Range& pr,
877 const Range& rr,
878 const Range& cr)
879 : ConstTensor4View(data, br, pr, rr, cr) {
880 // Nothing to do here.
881}
882
901Tensor4View::Tensor4View(Numeric* data,
902 const Range& pb,
903 const Range& pp,
904 const Range& pr,
905 const Range& pc,
906 const Range& nb,
907 const Range& np,
908 const Range& nr,
909 const Range& nc)
910 : ConstTensor4View(data, pb, pp, pr, pc, nb, np, nr, nc) {
911 // Nothing to do here.
912}
913
918void copy(ConstIterator4D origin,
919 const ConstIterator4D& end,
920 Iterator4D target) {
921 for (; origin != end; ++origin, ++target) {
922 // We use the copy function for the next smaller rank of tensor
923 // recursively:
924 copy(origin->begin(), origin->end(), target->begin());
925 }
926}
927
929void copy(Numeric x, Iterator4D target, const Iterator4D& end) {
930 for (; target != end; ++target) {
931 // We use the copy function for the next smaller rank of tensor
932 // recursively:
933 copy(x, target->begin(), target->end());
934 }
935}
936
937// Functions for Tensor4:
938// ---------------------
939
942Tensor4::Tensor4(Index b, Index p, Index r, Index c)
943 : Tensor4View(new Numeric[b * p * r * c],
944 Range(0, b, p * r * c),
945 Range(0, p, r * c),
946 Range(0, r, c),
947 Range(0, c)) {
948 // Nothing to do here.
949}
950
952Tensor4::Tensor4(Index b, Index p, Index r, Index c, Numeric fill)
953 : Tensor4View(new Numeric[b * p * r * c],
954 Range(0, b, p * r * c),
955 Range(0, p, r * c),
956 Range(0, r, c),
957 Range(0, c)) {
958 // Here we can access the raw memory directly, for slightly
959 // increased efficiency:
960 std::fill_n(mdata, b * p * r * c, fill);
961}
962
965Tensor4::Tensor4(const ConstTensor4View& m)
966 : Tensor4View(new Numeric[m.nbooks() * m.npages() * m.nrows() * m.ncols()],
967 Range(0, m.nbooks(), m.npages() * m.nrows() * m.ncols()),
968 Range(0, m.npages(), m.nrows() * m.ncols()),
969 Range(0, m.nrows(), m.ncols()),
970 Range(0, m.ncols())) {
971 copy(m.begin(), m.end(), begin());
972}
973
976Tensor4::Tensor4(const Tensor4& m)
977 : Tensor4View(new Numeric[m.nbooks() * m.npages() * m.nrows() * m.ncols()],
978 Range(0, m.nbooks(), m.npages() * m.nrows() * m.ncols()),
979 Range(0, m.npages(), m.nrows() * m.ncols()),
980 Range(0, m.nrows(), m.ncols()),
981 Range(0, m.ncols())) {
982 // There is a catch here: If m is an empty tensor, then it will have
983 // dimensions of size 0. But these are used to initialize the stride
984 // for higher dimensions! Thus, this method has to be consistent
985 // with the behaviour of Range::Range. For now, Range::Range allows
986 // also stride 0.
987 std::memcpy(mdata,
988 m.mdata,
989 nbooks() * npages() * nrows() * ncols() * sizeof(Numeric));
990}
991
993
1016Tensor4& Tensor4::operator=(const Tensor4& x) {
1017 if (this != &x) {
1018 resize(x.nbooks(), x.npages(), x.nrows(), x.ncols());
1019 std::memcpy(mdata,
1020 x.mdata,
1021 nbooks() * npages() * nrows() * ncols() * sizeof(Numeric));
1022 }
1023 return *this;
1024}
1025
1027Tensor4& Tensor4::operator=(Tensor4&& x) noexcept {
1028 if (this != &x) {
1029 delete[] mdata;
1030 mdata = x.mdata;
1031 mbr = x.mbr;
1032 mpr = x.mpr;
1033 mrr = x.mrr;
1034 mcr = x.mcr;
1035 x.mbr = Range(0, 0);
1036 x.mpr = Range(0, 0);
1037 x.mrr = Range(0, 0);
1038 x.mcr = Range(0, 0);
1039 x.mdata = nullptr;
1040 }
1041 return *this;
1042}
1043
1046Tensor4& Tensor4::operator=(Numeric x) {
1047 std::fill_n(mdata, nbooks() * npages() * nrows() * ncols(), x);
1048 return *this;
1049}
1050
1054void Tensor4::resize(Index b, Index p, Index r, Index c) {
1055 ARTS_ASSERT(0 <= b);
1056 ARTS_ASSERT(0 <= p);
1057 ARTS_ASSERT(0 <= r);
1058 ARTS_ASSERT(0 <= c);
1059
1060 if (mbr.mextent != b || mpr.mextent != p || mrr.mextent != r ||
1061 mcr.mextent != c) {
1062 delete[] mdata;
1063 mdata = new Numeric[b * p * r * c];
1064
1065 mbr.mstart = 0;
1066 mbr.mextent = b;
1067 mbr.mstride = p * r * c;
1068
1069 mpr.mstart = 0;
1070 mpr.mextent = p;
1071 mpr.mstride = r * c;
1072
1073 mrr.mstart = 0;
1074 mrr.mextent = r;
1075 mrr.mstride = c;
1076
1077 mcr.mstart = 0;
1078 mcr.mextent = c;
1079 mcr.mstride = 1;
1080 }
1081}
1082
1084void swap(Tensor4& t1, Tensor4& t2) {
1085 std::swap(t1.mbr, t2.mbr);
1086 std::swap(t1.mpr, t2.mpr);
1087 std::swap(t1.mrr, t2.mrr);
1088 std::swap(t1.mcr, t2.mcr);
1089 std::swap(t1.mdata, t2.mdata);
1090}
1091
1094Tensor4::~Tensor4() {
1095 // cout << "Destroying a Tensor4:\n"
1096 // << *this << "\n........................................\n";
1097 delete[] mdata;
1098}
1099
1115void transform(Tensor4View y, double (&my_func)(double), ConstTensor4View x) {
1116 // Check dimensions:
1117 ARTS_ASSERT(y.nbooks() == x.nbooks());
1118 ARTS_ASSERT(y.npages() == x.npages());
1119 ARTS_ASSERT(y.nrows() == x.nrows());
1120 ARTS_ASSERT(y.ncols() == x.ncols());
1121
1122 const ConstIterator4D xe = x.end();
1123 ConstIterator4D xi = x.begin();
1124 Iterator4D yi = y.begin();
1125 for (; xi != xe; ++xi, ++yi) {
1126 // Use the transform function of lower dimensional tensors
1127 // recursively:
1128 transform(*yi, my_func, *xi);
1129 }
1130}
1131
1133Numeric max(const ConstTensor4View& x) {
1134 const ConstIterator4D xe = x.end();
1135 ConstIterator4D xi = x.begin();
1136
1137 // Initial value for max:
1138 Numeric themax = max(*xi);
1139 ++xi;
1140
1141 for (; xi != xe; ++xi) {
1142 // Use the max function of lower dimensional tensors
1143 // recursively:
1144 Numeric maxi = max(*xi);
1145 if (maxi > themax) themax = maxi;
1146 }
1147
1148 return themax;
1149}
1150
1152Numeric min(const ConstTensor4View& x) {
1153 const ConstIterator4D xe = x.end();
1154 ConstIterator4D xi = x.begin();
1155
1156 // Initial value for min:
1157 Numeric themin = min(*xi);
1158 ++xi;
1159
1160 for (; xi != xe; ++xi) {
1161 // Use the min function of lower dimensional tensors
1162 // recursively:
1163 Numeric mini = min(*xi);
1164 if (mini < themin) themin = mini;
1165 }
1166
1167 return themin;
1168}
1169
1171// Helper function for debugging
1172#ifndef NDEBUG
1173
1189Numeric debug_tensor4view_get_elem(
1190 Tensor4View& tv, Index b, Index p, Index r, Index c) {
1191 return tv(b, p, r, c);
1192}
1193
1194#endif
ConstIterator4D
Const version of Iterator4D.
Definition: matpackIV.h:76
ConstIterator4D::operator*
const ConstTensor3View & operator*() const
Dereferencing.
Definition: matpackIV.cc:41
ConstIterator4D::operator->
const ConstTensor3View * operator->() const
The -> operator is needed, so that we can write i->begin() to get the 3D iterators.
Definition: matpackIV.cc:38
ConstIterator4D::msv
ConstTensor3View msv
Current position.
Definition: matpackIV.h:111
ConstMatrixView
A constant view of a Matrix.
Definition: matpackI.h:1050
ConstTensor3View
A constant view of a Tensor3.
Definition: matpackIII.h:130
ConstTensor3View::begin
ConstIterator3D begin() const
Return const iterator to first page.
Definition: matpackIII.cc:130
ConstTensor3View::end
ConstIterator3D end() const
Return const iterator behind last page.
Definition: matpackIII.cc:136
ConstTensor4View
A constant view of a Tensor4.
Definition: matpackIV.h:131
ConstTensor4View::ncols
Index ncols() const noexcept
Definition: matpackIV.h:142
ConstTensor4View::mrr
Range mrr
The row range of mdata that is actually used.
Definition: matpackIV.h:274
ConstTensor4View::ConstTensor4View
ConstTensor4View()=default
ConstTensor4View::mbr
Range mbr
The book range of mdata that is actually used.
Definition: matpackIV.h:270
ConstTensor4View::operator()
ConstTensor4View operator()(const Range &b, const Range &p, const Range &r, const Range &c) const
Const index operator for subrange.
Definition: matpackIV.cc:49
ConstTensor4View::size
Index size() const noexcept
Definition: matpackIV.h:145
ConstTensor4View::mcr
Range mcr
The column range of mdata that is actually used.
Definition: matpackIV.h:276
ConstTensor4View::mdata
Numeric * mdata
Pointer to the plain C array that holds the data.
Definition: matpackIV.h:278
ConstTensor4View::mpr
Range mpr
The page range of mdata that is actually used.
Definition: matpackIV.h:272
ConstTensor4View::end
ConstIterator4D end() const
Return const iterator behind last book.
Definition: matpackIV.cc:367
ConstTensor4View::begin
ConstIterator4D begin() const
Return const iterator to first book.
Definition: matpackIV.cc:361
ConstTensor4View::nrows
Index nrows() const noexcept
Definition: matpackIV.h:141
ConstTensor4View::nbooks
Index nbooks() const noexcept
Definition: matpackIV.h:139
ConstTensor4View::npages
Index npages() const noexcept
Definition: matpackIV.h:140
ConstVectorView
A constant view of a Vector.
Definition: matpackI.h:517
Iterator4D
Implementation of Tensors of Rank 4.
Definition: matpackIV.h:38
Iterator4D::operator->
Tensor3View * operator->()
The -> operator is needed, so that we can write i->begin() to get the 3D iterators.
Definition: matpackIV.cc:31
Iterator4D::operator*
Tensor3View & operator*()
Dereferencing.
Definition: matpackIV.cc:34
Iterator4D::msv
Tensor3View msv
Current position.
Definition: matpackIV.h:70
MatrixView
The MatrixView class.
Definition: matpackI.h:1169
Range
The range class.
Definition: matpackI.h:166
Range::mstart
Index mstart
The start index.
Definition: matpackI.h:368
Range::mstride
Index mstride
The stride.
Definition: matpackI.h:377
Range::mextent
Index mextent
The number of elements.
Definition: matpackI.h:375
Tensor3View
The Tensor3View class.
Definition: matpackIII.h:244
Tensor3View::begin
Iterator3D begin()
Return iterator to first page.
Definition: matpackIII.cc:329
Tensor3View::end
Iterator3D end()
Return iterator behind last page.
Definition: matpackIII.cc:334
Tensor4View
The Tensor4View class.
Definition: matpackIV.h:290
Tensor4View::begin
Iterator4D begin()
Return iterator to first book.
Definition: matpackIV.cc:696
Tensor4View::get_c_array
const Numeric * get_c_array() const ARTS_NOEXCEPT
Conversion to plain C-array.
Definition: matpackIV.cc:342
Tensor4View::operator+=
Tensor4View & operator+=(Numeric x)
Addition of scalar.
Definition: matpackIV.cc:783
Tensor4View::operator()
Tensor4View operator()(const Range &b, const Range &p, const Range &r, const Range &c)
Index operator for subrange.
Definition: matpackIV.cc:444
Tensor4View::operator-=
Tensor4View & operator-=(Numeric x)
Subtraction of scalar.
Definition: matpackIV.cc:792
Tensor4View::operator=
Tensor4View & operator=(const ConstTensor4View &v)
Assignment operator.
Definition: matpackIV.cc:713
Tensor4View::end
Iterator4D end()
Return iterator behind last book.
Definition: matpackIV.cc:702
Tensor4View::operator/=
Tensor4View & operator/=(Numeric x)
Division by scalar.
Definition: matpackIV.cc:774
Tensor4View::operator*=
Tensor4View & operator*=(Numeric x)
Multiplication by scalar.
Definition: matpackIV.cc:765
Tensor4View::Tensor4View
Tensor4View()=default
Tensor4
The Tensor4 class.
Definition: matpackIV.h:427
Tensor4::~Tensor4
virtual ~Tensor4()
Destructor for Tensor4.
Definition: matpackIV.cc:1094
Tensor4::resize
void resize(Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackIV.cc:1054
Tensor4::Tensor4
Tensor4()=default
Tensor4::operator=
Tensor4 & operator=(const Tensor4 &x)
Assignment operator from another tensor.
Definition: matpackIV.cc:1016
VectorView
The VectorView class.
Definition: matpackI.h:658
ARTS_NOEXCEPT
#define ARTS_NOEXCEPT
Definition: debug.h:80
ARTS_ASSERT
#define ARTS_ASSERT(condition,...)
Definition: debug.h:83
exceptions.h
The declarations of all the exception classes.
eb
#define eb
Definition: legacy_continua.cc:22126
copy
void copy(ConstIterator4D origin, const ConstIterator4D &end, Iterator4D target)
Copy data between begin and end to target.
Definition: matpackIV.cc:918
debug_tensor4view_get_elem
Numeric debug_tensor4view_get_elem(Tensor4View &tv, Index b, Index p, Index r, Index c)
Helper function to access tensor elements.
Definition: matpackIV.cc:1189
max
Numeric max(const ConstTensor4View &x)
Max function, tensor version.
Definition: matpackIV.cc:1133
min
Numeric min(const ConstTensor4View &x)
Min function, tensor version.
Definition: matpackIV.cc:1152
operator<<
std::ostream & operator<<(std::ostream &os, const ConstTensor4View &v)
Output operator.
Definition: matpackIV.cc:420
transform
void transform(Tensor4View y, double(&my_func)(double), ConstTensor4View x)
A generic transform function for tensors, which can be used to implement mathematical functions opera...
Definition: matpackIV.cc:1115
swap
void swap(Tensor4 &t1, Tensor4 &t2)
Swaps two objects.
Definition: matpackIV.cc:1084
Numeric
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
Index
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
Zeeman::end
constexpr Rational end(Rational Ju, Rational Jl, Polarization type) noexcept
Gives the largest M for a polarization type of this transition.
Definition: zeemandata.h:113
v
#define v
a
#define a
c
#define c
b
#define b