Changeset 1757 for branches/alma/src/MathUtils.cpp
- Timestamp:
- 06/09/10 19:03:06 (14 years ago)
- Location:
- branches/alma
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/alma
-
Property
svn:ignore
set to
.sconf_temp
.sconsign.dblite
-
Property
svn:mergeinfo
set to
/branches/asap-3.x merged eligible
-
Property
svn:ignore
set to
-
branches/alma/src/MathUtils.cpp
r1603 r1757 39 39 #include <scimath/Mathematics/MedianSlider.h> 40 40 #include <casa/Exceptions/Error.h> 41 42 #include <scimath/Fitting/LinearFit.h> 43 #include <scimath/Functionals/Polynomial.h> 44 #include <scimath/Mathematics/AutoDiff.h> 45 41 46 42 47 #include "MathUtils.h" … … 183 188 MedianSlider ms(hwidth); 184 189 Slice sl(0, fwidth-1); 185 Float medval = ms.add(const_cast<Vector<Float>& >(in)(sl), 190 Float medval = ms.add(const_cast<Vector<Float>& >(in)(sl), 186 191 const_cast<Vector<Bool>& >(flag)(sl)); 187 192 uInt n = in.nelements(); 188 193 for (uInt i=hwidth; i<(n-hwidth); ++i) { 189 194 // add data value 190 out[i] = ms.add(in[i+hwidth], flag[i+hwidth]); 191 outflag[i] = (ms.nval() == 0); 192 } 193 // replicate edge values from fi srt value with full width of values195 out[i] = ms.add(in[i+hwidth], flag[i+hwidth]); 196 outflag[i] = (ms.nval() == 0); 197 } 198 // replicate edge values from first value with full width of values 194 199 for (uInt i=0;i<hwidth;++i) { 195 200 out[i] = out[hwidth]; 196 outflag[i] = outflag[hwidth]; 201 outflag[i] = outflag[hwidth]; 197 202 out[n-1-i] = out[n-1-hwidth]; 198 outflag[n-1-i] = outflag[n-1-hwidth]; 199 } 200 } 203 outflag[n-1-i] = outflag[n-1-hwidth]; 204 } 205 } 206 207 void mathutil::polyfit(Vector<Float>& out, Vector<Bool>& outmask, 208 const Vector<Float>& in, const Vector<Bool>& mask, 209 float width, int order) 210 { 211 Int hwidth = Int(width+0.5); 212 Int fwidth = hwidth*2+1; 213 out.resize(in.nelements()); 214 outmask.resize(mask.nelements()); 215 LinearFit<Float> fitter; 216 Polynomial<Float> poly(order); 217 fitter.setFunction(poly); 218 Vector<Float> sigma(fwidth); 219 sigma = 1.0; 220 Vector<Float> parms; 221 Vector<Float> x(fwidth); 222 indgen(x); 223 224 uInt n = in.nelements(); 225 226 for (uInt i=hwidth; i<(n-hwidth); ++i) { 227 // add data value 228 if (mask[i]) { 229 Slice sl(i-hwidth, fwidth); 230 const Vector<Float> &y = const_cast<Vector<Float>& >(in)(sl); 231 const Vector<Bool> &m = const_cast<Vector<Bool>& >(mask)(sl); 232 parms = fitter.fit(x, y, sigma, &m); 233 234 poly.setCoefficients(parms); 235 out[i] = poly(x[hwidth]);//cout << in[i] <<"->"<<out[i]<<endl; 236 } else { 237 out[i] = in[i]; 238 } 239 outmask[i] = mask[i]; 240 } 241 // replicate edge values from first value with full width of values 242 for (uInt i=0;i<hwidth;++i) { 243 out[i] = out[hwidth]; 244 outmask[i] = outmask[hwidth]; 245 out[n-1-i] = out[n-1-hwidth]; 246 outmask[n-1-i] = outmask[n-1-hwidth]; 247 } 248 }
Note: See TracChangeset
for help on using the changeset viewer.