source: branches/Release12/src/MathUtils2.cc@ 1051

Last change on this file since 1051 was 452, checked in by mar637, 20 years ago
  • Added extendLastArrayAxis
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 4.0 KB
Line 
1//#---------------------------------------------------------------------------
2//# MathUtilities.cc: General math operations
3//#---------------------------------------------------------------------------
4//# Copyright (C) 2004
5//# ATNF
6//#
7//# This program is free software; you can redistribute it and/or modify it
8//# under the terms of the GNU General Public License as published by the Free
9//# Software Foundation; either version 2 of the License, or (at your option)
10//# any later version.
11//#
12//# This program is distributed in the hope that it will be useful, but
13//# WITHOUT ANY WARRANTY; without even the implied warranty of
14//# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
15//# Public License for more details.
16//#
17//# You should have received a copy of the GNU General Public License along
18//# with this program; if not, write to the Free Software Foundation, Inc.,
19//# 675 Massachusetts Ave, Cambridge, MA 02139, USA.
20//#
21//# Correspondence concerning this software should be addressed as follows:
22//# Internet email: Malte.Marquarding@csiro.au
23//# Postal address: Malte Marquarding,
24//# Australia Telescope National Facility,
25//# P.O. Box 76,
26//# Epping, NSW, 2121,
27//# AUSTRALIA
28//#
29//# $Id:
30//#---------------------------------------------------------------------------
31
32#include <casa/aips.h>
33#include <casa/Arrays/Array.h>
34#include <casa/Arrays/ArrayIter.h>
35#include <casa/Arrays/Vector.h>
36#include <casa/Arrays/VectorSTLIterator.h>
37
38#include "MathUtils.h"
39
40using namespace casa;
41
42template <class T>
43void mathutil::hanning(Vector<T>& out, Vector<Bool>& outmask,
44 const Vector<T>& in, const Vector<Bool>& mask,
45 Bool relaxed, Bool ignoreOther) {
46
47 Vector< Vector<T> > weights(8);
48 Vector<Float> vals(3);
49 vals = 0.0;weights[0] = vals;// FFF
50 vals[0] = 1.0; vals[1] = 0.0; vals[2] = 0.0; weights[1] = vals;// TFF
51 vals[0] = 0.0; vals[1] = 1.0; vals[2] = 0.0; weights[2] = vals;// FTF
52 vals[0] = 1.0/3.0; vals[1] = 2.0/3.0; vals[2] = 0.0; weights[3] = vals;// TTF
53 vals[0] = 0.0; vals[1] = 0.0; vals[2] = 1.0;weights[4] = vals;// FFT
54 vals[0] = 0.5; vals[1] = 0.0; vals[2] = 0.5; weights[5] = vals;// TFT
55 vals[0] = 0.0; vals[1] = 2.0/3.0; vals[2] = 1.0/3.0; weights[6] = vals;// FTT
56 vals[0] = 0.25; vals[1] = 0.5; vals[2] = 0.25; weights[7] = vals;// TTT
57 // Chris' case
58 Vector<Bool> weighted(8);
59 if (relaxed) {
60 weighted = False;
61 weighted[7] = True;
62
63 } else {
64 weighted = True;
65 weighted[0] = False;
66 }
67
68 out.resize(in.nelements());
69 outmask.resize(mask.nelements());
70
71 // make special case for first and last
72 /// ...here
73 // loop from 1..n-2
74 uInt i = 1;
75 VectorSTLIterator<T> outit(out);
76 outit++;
77 VectorSTLIterator<Bool> outmit(outmask);outmit++;
78 uInt m;Vector<T>* w;
79 for (VectorSTLIterator<T> it = outit;it != out.end()-1;++it) {
80
81 m = mask[i-1] + 2*mask[i] + 4*mask[i+1];
82 w = &(weights[m]);
83 if (weighted[m]) {
84 (*it) = (*w)[0]*in[i-1] + (*w)[1]*in[i] + (*w)[2]*in[i+1];
85 (*outmit) = True;
86 } else { // mask it
87 (*outmit) = False;
88 }
89 ++i;
90 ++outmit;
91 }
92}
93
94template <class T>
95uInt mathutil::addEntry (Vector<T>& list, T val)
96{
97 uInt n = list.nelements();
98 if (n>0) {
99 for (uInt i=0; i<n; i++) {
100 if (near(list[i],val)) {
101 return i;
102 }
103 }
104 }
105//
106 n++;
107 list.resize(n,True);
108 list(n-1) = val;
109 return n-1;
110}
111
112template <class T>
113void mathutil::extendLastArrayAxis(casa::Array<T>& out,
114 const casa::Array<T>& in,
115 const T& initVal)
116{
117
118 IPosition ipin = in.shape();
119 IPosition ipout = ipin; // copy?
120 // extend the axis by 1
121 uInt axis = in.ndim()-1;
122 ipout[axis] = ipin[axis]+1;
123 out.resize(ipout);
124 out = initVal;
125 ArrayIterator<T> itout(out,axis);
126 ReadOnlyArrayIterator<T> itin(in,axis);
127 while ( !itin.pastEnd() ) {
128 itout.array() = itin.array(); // copy vector by vector
129 itin.next(); itout.next();
130 }
131}
Note: See TracBrowser for help on using the repository browser.