Changeset 1924
- Timestamp:
- 09/14/10 12:17:16 (14 years ago)
- Location:
- branches/polybatch
- Files:
-
- 11 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/polybatch/python/scantable.py
r1920 r1924 1158 1158 return msk 1159 1159 1160 def get_masklist(self, mask=None, row=0 ):1160 def get_masklist(self, mask=None, row=0, silent=False): 1161 1161 """\ 1162 1162 Compute and return a list of mask windows, [min, max]. … … 1192 1192 if not i: 1193 1193 msg += "\nThis mask is only valid for IF=%d" % (self.getif(i)) 1194 asaplog.push(msg) 1194 if not silent: 1195 asaplog.push(msg) 1195 1196 masklist=[] 1196 1197 ist, ien = None, None … … 1934 1935 workscan._setspectrum(f.fitter.getresidual(), r) 1935 1936 self.blpars.append(f.get_parameters()) 1936 self.masklists.append(workscan.get_masklist(f.mask, row=r ))1937 self.masklists.append(workscan.get_masklist(f.mask, row=r, silent=True)) 1937 1938 self.actualmask.append(f.mask) 1938 1939 … … 1949 1950 raise RuntimeError(msg) 1950 1951 1951 1952 def poly_baseline(self, mask=None, order=0, plot=False, batch=False, insitu=None, rows=None): 1952 @asaplog_post_dec 1953 def poly_baseline(self, mask=None, order=0, plot=False, batch=False, 1954 insitu=None, rows=None): 1953 1955 """\ 1954 1956 Return a scan which has been baselined (all rows) by a polynomial. … … 1973 1975 """ 1974 1976 if insitu is None: insitu = rcParams["insitu"] 1977 varlist = vars() 1975 1978 if insitu: 1976 1979 workscan = self … … 1978 1981 workscan = self.copy() 1979 1982 1980 varlist = vars()1981 1983 nchan = workscan.nchan() 1982 1984 … … 1991 1993 1992 1994 if len(rows) > 0: 1993 self.blpars = []1994 self.masklists = []1995 self.actualmask = []1995 workscan.blpars = [] 1996 workscan.masklists = [] 1997 workscan.actualmask = [] 1996 1998 1997 1999 if batch: 1998 for r in rows: 1999 workscan._poly_baseline_batch(mask, order, r) 2000 workscan._poly_baseline_batch(mask, order) 2000 2001 elif plot: 2001 2002 f = fitter() … … 2016 2017 continue 2017 2018 workscan._setspectrum(f.fitter.getresidual(), r) 2018 self.blpars.append(f.get_parameters())2019 self.masklists.append(workscan.get_masklist(f.mask, row=r))2020 self.actualmask.append(f.mask)2019 workscan.blpars.append(f.get_parameters()) 2020 workscan.masklists.append(workscan.get_masklist(f.mask, row=r)) 2021 workscan.actualmask.append(f.mask) 2021 2022 2022 2023 f._p.unmap() 2023 2024 f._p = None 2024 2025 else: 2025 import array2026 2026 for r in rows: 2027 pars = array.array("f", [0.0 for i in range(order+1)]) 2028 pars_adr = pars.buffer_info()[0] 2029 pars_len = pars.buffer_info()[1] 2027 fitparams = workscan._poly_baseline(mask, order, r) 2028 params = fitparams.getparameters() 2029 fmtd = ", ".join(["p%d = %3.6f" % (i, v) for i, v in enumerate(params)]) 2030 errors = fitparams.geterrors() 2031 fmask = mask_and(mask, workscan._getmask(r)) 2032 2033 workscan.blpars.append({"params":params, 2034 "fixed": fitparams.getfixedparameters(), 2035 "formatted":fmtd, "errors":errors}) 2036 workscan.masklists.append(workscan.get_masklist(fmask, r, silent=True)) 2037 workscan.actualmask.append(fmask) 2030 2038 2031 errs = array.array("f", [0.0 for i in range(order+1)]) 2032 errs_adr = errs.buffer_info()[0] 2033 errs_len = errs.buffer_info()[1] 2034 2035 fmsk = array.array("i", [1 for i in range(nchan)]) 2036 fmsk_adr = fmsk.buffer_info()[0] 2037 fmsk_len = fmsk.buffer_info()[1] 2038 2039 workscan._poly_baseline(mask, order, r, pars_adr, pars_len, errs_adr, errs_len, fmsk_adr, fmsk_len) 2040 2041 params = pars.tolist() 2042 fmtd = "" 2043 for i in xrange(len(params)): fmtd += " p%d= %3.6f," % (i, params[i]) 2044 fmtd = fmtd[:-1] # remove trailing "," 2045 errors = errs.tolist() 2046 fmask = fmsk.tolist() 2047 for i in xrange(len(fmask)): fmask[i] = (fmask[i] > 0) # transform (1/0) -> (True/False) 2048 2049 self.blpars.append({"params":params, "fixed":[], "formatted":fmtd, "errors":errors}) 2050 self.masklists.append(workscan.get_masklist(fmask, r)) 2051 self.actualmask.append(fmask) 2052 2053 asaplog.push(str(fmtd)) 2039 asaplog.push(fmtd) 2054 2040 2055 2041 workscan._add_history("poly_baseline", varlist) … … 2202 2188 2203 2189 # Show mask list 2204 masklist=workscan.get_masklist(f.mask, row=r )2190 masklist=workscan.get_masklist(f.mask, row=r, silent=True) 2205 2191 msg = "mask range: "+str(masklist) 2206 2192 asaplog.push(msg, False) -
branches/polybatch/src/STFit.cpp
r1000 r1924 70 70 table_.addColumn(ArrayColumnDesc<Int>("COMPONENTS")); 71 71 table_.addColumn(ArrayColumnDesc<Double>("PARAMETERS")); 72 // table_.addColumn(ArrayColumnDesc<Double>("ERRORS")); 72 73 table_.addColumn(ArrayColumnDesc<Bool>("PARMASKS")); 73 74 table_.addColumn(ArrayColumnDesc<String>("FRAMEINFO")); … … 77 78 compCol_.attach(table_,"COMPONENTS"); 78 79 parCol_.attach(table_,"PARAMETERS"); 80 // errCol_.attach(table_,"ERRORS"); 79 81 maskCol_.attach(table_,"PARMASKS"); 80 82 frameCol_.attach(table_,"FRAMEINFO"); … … 102 104 // add new row if new id 103 105 if ( !foundentry ) table_.addRow(); 106 104 107 funcCol_.put(rno, mathutil::toVectorString(fit.getFunctions())); 105 108 compCol_.put(rno, Vector<Int>(fit.getComponents())); 106 parCol_.put(rno, Vector<Double>(fit.getParameters())); 109 const std::vector<float>& pvec = fit.getParameters(); 110 Vector<Double> dvec(pvec.size()); 111 for (size_t i=0; i < dvec.nelements(); ++i) { 112 dvec[i] = Double(pvec[i]); 113 } 114 parCol_.put(rno, dvec); 115 /* 116 const std::vector<double>& evec = fit.getErrors(); 117 for (size_t i=0; i < dvec.nelements(); ++i) { 118 dvec[i] = Double(evec[i]); 119 } 120 errCol_.put(rno, dvec); 121 */ 107 122 maskCol_.put(rno, Vector<Bool>(fit.getParmasks())); 108 123 frameCol_.put(rno, mathutil::toVectorString(fit.getFrameinfo())); 109 124 idCol_.put(rno, resultid); 125 110 126 return resultid; 111 127 } … … 130 146 fit.setComponents(istl); 131 147 Vector<Double> dvec; 132 std::vector<double> dstl;133 148 rec.get("PARAMETERS", dvec); 149 std::vector<float> dstl(dvec.begin(), dvec.end()); 150 fit.setParameters(dstl); 151 /* 134 152 dvec.tovector(dstl); 135 153 fit.setParameters(dstl); 154 dvec.resize(); 155 rec.get("ERRORS", dvec); 156 dvec.tovector(dstl); 157 fit.setErrors(dstl); 158 */ 136 159 Vector<Bool> bvec; 137 160 std::vector<bool> bstl; -
branches/polybatch/src/STFit.h
r1353 r1924 49 49 casa::ArrayColumn<casa::Int> compCol_; 50 50 casa::ArrayColumn<casa::Double> parCol_; 51 // casa::ArrayColumn<casa::Double> errCol_; 51 52 casa::ArrayColumn<casa::Bool> maskCol_; 52 53 casa::ArrayColumn<casa::String> frameCol_; -
branches/polybatch/src/STFitEntry.cpp
r972 r1924 11 11 // 12 12 #include "STFitEntry.h" 13 #include <casa/iostream.h> 13 14 15 using namespace casa; 14 16 namespace asap { 15 17 16 STFitEntry::STFitEntry() 18 STFitEntry::STFitEntry() 19 /* : 20 functions_( std::vector<std::string>()), 21 components_(std::vector<int>()), 22 parameters_(std::vector<float>()), 23 errors_(std::vector<float>()), 24 parmasks_(std::vector<bool>()), 25 frameinfo_(std::vector<std::string>()) 26 */ 17 27 { 18 28 } … … 20 30 { 21 31 if ( this != &other ) { 22 this->functions_ = other.functions_;32 this->functions_ = std::vector<std::string>(); 23 33 this->components_ = other.components_; 24 34 this->parameters_ = other.parameters_; 25 this-> parmasks_ = other.parmasks_;35 this->errors_ = other.errors_; 26 36 this->frameinfo_ = other.frameinfo_; 27 37 } 28 38 } 29 39 40 STFitEntry& STFitEntry::operator=(const STFitEntry& other) 41 { 42 if ( this != &other ) { 43 this->functions_ = other.functions_; 44 this->components_ = other.components_; 45 this->parameters_ = other.parameters_; 46 this->errors_ = other.errors_; 47 this->parmasks_ = other.parmasks_; 48 this->frameinfo_ = other.frameinfo_; 49 } 50 return *this; 51 } 30 52 31 53 STFitEntry::~STFitEntry() -
branches/polybatch/src/STFitEntry.h
r972 r1924 28 28 STFitEntry(const STFitEntry& other); 29 29 30 STFitEntry& operator=(const STFitEntry& other); 31 30 32 ~STFitEntry(); 31 33 … … 34 36 void setComponents(const std::vector<int>& c) 35 37 { components_ = c; } 36 void setParameters(const std::vector< double>& p)38 void setParameters(const std::vector<float>& p) 37 39 { parameters_ = p; } 40 void setErrors(const std::vector<float>& p) 41 { errors_ = p; } 38 42 void setParmasks(const std::vector<bool>& m) 39 43 { parmasks_ = m; } 40 44 void setFrameinfo(const std::vector<std::string>& f) 41 45 { frameinfo_ = f; } 42 43 46 std::vector<std::string> getFunctions() const { return functions_; } 44 47 std::vector<int> getComponents() const { return components_; } 45 std::vector<double> getParameters() const { return parameters_; } 48 std::vector<float> getParameters() const { return parameters_; } 49 std::vector<float> getErrors() const { return errors_; } 46 50 std::vector<bool> getParmasks() const { return parmasks_; } 47 51 std::vector<std::string> getFrameinfo() const { return frameinfo_; } 48 52 49 53 private: 54 50 55 std::vector<std::string> functions_; 51 56 std::vector<int> components_; 52 std::vector<double> parameters_; 57 std::vector<float> parameters_; 58 std::vector<float> errors_; 53 59 std::vector<bool> parmasks_; 54 60 std::vector<std::string> frameinfo_; -
branches/polybatch/src/STFitter.cpp
r1819 r1924 142 142 if (ncomp < 1) throw (AipsError("Need at least one gaussian to fit.")); 143 143 funcs_.resize(ncomp); 144 funcnames_.clear(); 145 funccomponents_.clear(); 144 146 for (Int k=0; k<ncomp; ++k) { 145 147 funcs_[k] = new Gaussian1D<Float>(); 148 funcnames_.push_back(expr); 149 funccomponents_.push_back(3); 146 150 } 147 151 } else if (expr == "poly") { 148 152 funcs_.resize(1); 153 funcnames_.clear(); 154 funccomponents_.clear(); 149 155 funcs_[0] = new Polynomial<Float>(ncomp); 156 funcnames_.push_back(expr); 157 funccomponents_.push_back(ncomp); 150 158 } else if (expr == "lorentz") { 151 159 if (ncomp < 1) throw (AipsError("Need at least one lorentzian to fit.")); 152 160 funcs_.resize(ncomp); 161 funcnames_.clear(); 162 funccomponents_.clear(); 153 163 for (Int k=0; k<ncomp; ++k) { 154 164 funcs_[k] = new Lorentzian1D<Float>(); 165 funcnames_.push_back(expr); 166 funccomponents_.push_back(3); 155 167 } 156 168 } else { … … 409 421 } 410 422 423 STFitEntry Fitter::getFitEntry() const 424 { 425 STFitEntry fit; 426 fit.setParameters(getParameters()); 427 fit.setErrors(getErrors()); 428 fit.setComponents(funccomponents_); 429 fit.setFunctions(funcnames_); 430 fit.setParmasks(getFixedParameters()); 431 return fit; 432 } -
branches/polybatch/src/STFitter.h
r1391 r1924 40 40 #include <scimath/Functionals/CompoundFunction.h> 41 41 42 #include "STFitEntry.h" 43 42 44 namespace asap { 43 45 … … 69 71 70 72 std::vector<float> evaluate(int whichComp) const; 73 74 STFitEntry getFitEntry() const; 75 71 76 private: 72 77 void clear(); … … 75 80 casa::Vector<casa::Bool> m_; 76 81 casa::PtrBlock<casa::Function<casa::Float>* > funcs_; 82 std::vector<std::string> funcnames_; 83 std::vector<int> funccomponents_; 84 77 85 //Bool estimateSet_; 78 86 casa::Float chisquared_; -
branches/polybatch/src/Scantable.cpp
r1919 r1924 1750 1750 } 1751 1751 1752 void Scantable::polyBaselineBatch(const std::vector<bool>& mask, int order, int rowno) 1752 void Scantable::polyBaselineBatch(const std::vector<bool>& mask, int order) 1753 { 1754 Fitter fitter = Fitter(); 1755 for (uInt rowno=0; rowno < nrow(); ++rowno) { 1756 doPolyBaseline(mask, order, rowno, fitter); 1757 setSpectrum(fitter.getResidual(), rowno); 1758 } 1759 } 1760 1761 STFitEntry Scantable::polyBaseline(const std::vector<bool>& mask, int order, int rowno) 1753 1762 { 1754 1763 Fitter fitter = Fitter(); 1755 1764 doPolyBaseline(mask, order, rowno, fitter); 1756 1765 setSpectrum(fitter.getResidual(), rowno); 1757 } 1758 1759 void Scantable::polyBaseline(const std::vector<bool>& mask, int order, int rowno, long pars_ptr, long pars_size, long errs_ptr, long errs_size, long fmask_ptr, long fmask_size) 1760 { 1761 Fitter fitter = Fitter(); 1762 doPolyBaseline(mask, order, rowno, fitter); 1763 setSpectrum(fitter.getResidual(), rowno); 1764 1765 std::vector<float> pars = fitter.getParameters(); 1766 if (pars_size != pars.size()) { 1767 throw(AipsError("wrong pars size")); 1768 } 1769 float *ppars = reinterpret_cast<float*>(pars_ptr); 1770 for (int i = 0; i < pars_size; i++) { 1771 ppars[i] = pars[i]; 1772 } 1773 1774 std::vector<float> errs = fitter.getErrors(); 1775 if (errs_size != errs.size()) { 1776 throw(AipsError("wrong errors size")); 1777 } 1778 float *perrs = reinterpret_cast<float*>(errs_ptr); 1779 for (int i = 0; i < errs_size; i++) { 1780 perrs[i] = errs[i]; 1781 } 1782 1783 std::vector<bool> fmask = getMask(rowno); 1784 if (fmask_size != fmask.size()) { 1785 throw(AipsError("wrong fmask size")); 1786 } 1787 int *pfmask = reinterpret_cast<int*>(fmask_ptr); 1788 for (int i = 0; i < fmask_size; i++) { 1789 pfmask[i] = ((fmask[i] && mask[i]) ? 1 : 0); 1790 } 1766 return fitter.getFitEntry(); 1791 1767 } 1792 1768 -
branches/polybatch/src/Scantable.h
r1919 r1924 18 18 // AIPS++ 19 19 #include <casa/aips.h> 20 #include <casa/Containers/Record.h> 20 21 #include <casa/Arrays/MaskedArray.h> 21 22 #include <casa/BasicSL/String.h> … … 489 490 bool getFlagtraFast(int whichrow); 490 491 491 void polyBaselineBatch(const std::vector<bool>& mask, int order, int rowno); 492 void polyBaseline(const std::vector<bool>& mask, int order, int rowno, long pars_ptr, long pars_size, long errs_ptr, long errs_size, long fmask_ptr, long fmask_size); 493 492 void polyBaselineBatch(const std::vector<bool>& mask, int order); 493 STFitEntry polyBaseline(const std::vector<bool>& mask, int order, int rowno); 494 494 495 495 private: -
branches/polybatch/src/ScantableWrapper.h
r1919 r1924 19 19 #include "MathUtils.h" 20 20 #include "STFit.h" 21 #include "STFitEntry.h" 21 22 #include "Scantable.h" 22 23 #include "STCoordinate.h" … … 250 251 { table_->reshapeSpectrum( nmin, nmax ); } 251 252 252 void polyBaseline(const std::vector<bool>& mask, int order, int rowno, long pars_ptr, long pars_size, long errs_ptr, long errs_size, long fmask_ptr, long fmask_size)253 { table_->polyBaseline(mask, order, rowno, pars_ptr, pars_size, errs_ptr, errs_size, fmask_ptr, fmask_size); }254 255 void polyBaselineBatch(const std::vector<bool>& mask, int order , int rowno)256 { table_->polyBaselineBatch(mask, order , rowno); }253 STFitEntry polyBaseline(const std::vector<bool>& mask, int order, int rowno) 254 { return table_->polyBaseline(mask, order, rowno); } 255 256 void polyBaselineBatch(const std::vector<bool>& mask, int order) 257 { table_->polyBaselineBatch(mask, order); } 257 258 258 259 bool getFlagtraFast(int whichrow=0) const -
branches/polybatch/src/python_STFitEntry.cpp
r972 r1924 44 44 .def("getfixedparameters", &STFitEntry::getParmasks) 45 45 .def("getparameters", &STFitEntry::getParameters) 46 .def("geterrors", &STFitEntry::getErrors) 46 47 .def("getfunctions", &STFitEntry::getFunctions) 47 48 .def("getcomponents", &STFitEntry::getComponents) … … 49 50 .def("setfixedparameters", &STFitEntry::setParmasks) 50 51 .def("setparameters", &STFitEntry::setParameters) 52 .def("seterrors", &STFitEntry::setErrors) 51 53 .def("setfunctions", &STFitEntry::setFunctions) 52 54 .def("setcomponents", &STFitEntry::setComponents)
Note:
See TracChangeset
for help on using the changeset viewer.