source: trunk/src/SDMath.cc@ 23

Last change on this file since 23 was 15, checked in by mmarquar, 20 years ago

Added multiply function.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 5.5 KB
Line 
1//#---------------------------------------------------------------------------
2//# SDMath.cc: A collection of single dish mathematical operations
3//#---------------------------------------------------------------------------
4//# Copyright (C) 2004
5//# Malte Marquarding, 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#include <aips/aips.h>
32#include <aips/Utilities/String.h>
33#include <aips/Arrays/IPosition.h>
34#include <aips/Arrays/Array.h>
35#include <aips/Arrays/ArrayMath.h>
36#include <aips/Arrays/ArrayLogical.h>
37#include <aips/Arrays/MaskedArray.h>
38#include <aips/Arrays/MaskArrMath.h>
39#include <aips/Arrays/MaskArrLogi.h>
40
41#include <aips/Tables/Table.h>
42#include <aips/Tables/ScalarColumn.h>
43#include <aips/Tables/ArrayColumn.h>
44
45#include "SDContainer.h"
46#include "SDMemTable.h"
47
48#include "SDMath.h"
49
50using namespace atnf_sd;
51
52static CountedPtr<SDMemTable> SDMath::average(const CountedPtr<SDMemTable>& in) {
53 Table t = in->table();
54 ROArrayColumn<Float> tsys(t, "TSYS");
55 ROScalarColumn<Double> mjd(t, "TIME");
56 ROScalarColumn<String> srcn(t, "SRCNAME");
57 ROScalarColumn<Double> integr(t, "INTERVAL");
58 IPosition ip = in->rowAsMaskedArray(0).shape();
59 Array<Float> outarr(ip); outarr =0.0;
60 Array<Float> narr(ip);narr = 0.0;
61 Array<Float> narrinc(ip);narrinc = 1.0;
62
63 Array<Float> tsarr(tsys.shape(0));
64 Array<Float> outtsarr(tsys.shape(0));
65 Double tme = 0.0;
66 Double inttime = 0.0;
67
68 for (uInt i=0; i < t.nrow(); i++) {
69 // data stuff
70 MaskedArray<Float> marr(in->rowAsMaskedArray(i));
71 outarr += marr;
72 MaskedArray<Float> n(narrinc,marr.getMask());
73 narr += n;
74 // get
75 tsys.get(i, tsarr);// this is probably unneccessary as tsys should
76 outtsarr += tsarr; // be constant
77 Double tmp;
78 mjd.get(i,tmp);
79 tme += tmp;// average time
80 integr.get(i,tmp);
81 inttime += tmp;
82 }
83 // averaging using mask
84 MaskedArray<Float> nma(narr,(narr > Float(0)));
85 outarr /= nma;
86 Array<Bool> outflagsb = !(nma.getMask());
87 Array<uChar> outflags(outflagsb.shape());
88 convertArray(outflags,outflagsb);
89 SDContainer sc(ip(0),ip(1),ip(2),ip(3));
90
91 Int n = t.nrow();
92 outtsarr /= Float(n/2);
93 sc.timestamp = tme/Double(n/2);
94 sc.interval =inttime;
95 String tstr; srcn.getScalar(0,tstr);// get sourcename of "mid" point
96 sc.sourcename = tstr;
97 sc.putSpectrum(outarr);
98 sc.putFlags(outflags);
99 SDMemTable* sdmt = new SDMemTable(*in,True);
100 sdmt->putSDContainer(sc);
101 return CountedPtr<SDMemTable>(sdmt);
102}
103
104static CountedPtr<SDMemTable>
105SDMath::quotient(const CountedPtr<SDMemTable>& on,
106 const CountedPtr<SDMemTable>& off) {
107
108 Table ton = on->table();
109 Table toff = off->table();
110 ROArrayColumn<Float> tsys(toff, "TSYS");
111 ROScalarColumn<Double> mjd(ton, "TIME");
112 ROScalarColumn<Double> integr(ton, "INTERVAL");
113 ROScalarColumn<String> srcn(ton, "SRCNAME");
114 MaskedArray<Float> mon(on->rowAsMaskedArray(0));
115 MaskedArray<Float> moff(off->rowAsMaskedArray(0));
116 IPosition ipon = mon.shape();
117 IPosition ipoff = moff.shape();
118 Array<Float> tsarr(tsys.shape(0));
119
120 if (ipon != ipoff && ipon != tsarr.shape())
121 cerr << "on/off not conformant" << endl;
122
123 //IPosition test = mon.shape()/2;
124 MaskedArray<Float> tmp = (mon-moff);
125 Array<Float> out(tmp.getArray());
126 out /= moff;
127 out *= tsarr;
128 Array<Bool> outflagsb = !(mon.getMask() && moff.getMask());
129 Array<uChar> outflags(outflagsb.shape());
130 convertArray(outflags,outflagsb);
131
132 SDContainer sc(ipon(0),ipon(1),ipon(2),ipon(3));
133 String tstr; srcn.getScalar(0,tstr);// get sourcename of "on" scan
134 sc.sourcename = tstr;
135 Double tme; mjd.getScalar(0,tme);// get time of "on" scan
136 sc.timestamp = tme;
137 integr.getScalar(0,tme);
138 sc.interval = tme;
139 sc.putSpectrum(out);
140 sc.putFlags(outflags);
141 SDMemTable* sdmt = new SDMemTable(*on, True);
142 sdmt->putSDContainer(sc);
143 return CountedPtr<SDMemTable>(sdmt);
144
145}
146static CountedPtr<SDMemTable>
147SDMath::multiply(const CountedPtr<SDMemTable>& in, Float factor) {
148 SDMemTable* sdmt = new SDMemTable(*in);
149 Table t = sdmt->table();
150 ArrayColumn<Float> spec(t,"SPECTRA");
151
152 for (uInt i=0; i < t.nrow(); i++) {
153 // data stuff
154 MaskedArray<Float> marr(sdmt->rowAsMaskedArray(i));
155 marr *= factor;
156 spec.put(i, marr.getArray());
157 }
158 return CountedPtr<SDMemTable>(sdmt);
159}
160/*
161static Float SDMath::rms(const SDMemTable& in, uInt whichRow) {
162 Table t = in.table();
163 MaskedArray<Float> marr(in.rowAsMaskedArray(whichRow,True));
164
165}
166*/
Note: See TracBrowser for help on using the repository browser.