- Timestamp:
- 03/08/06 11:57:12 (19 years ago)
- Location:
- trunk/python
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/python/__init__.py
r876 r880 253 253 254 254 from asapmath import * 255 #from scantable import *256 #from asaplinefind import *255 from scantable import * 256 from asaplinefind import * 257 257 #from asapfit import * 258 258 -
trunk/python/asapfitter.py
r876 r880 452 452 else: 453 453 scan = self.data 454 rows = range(scan.nrow())454 rows = xrange(scan.nrow()) 455 455 from asap import asaplog 456 456 asaplog.push("Fitting:") 457 457 for r in rows: 458 out = " Scan[%d] Beam[%d] IF[%d] Pol[%d] Cycle[%d]" % (self.data.getscan(r),self.data.getbeam(r),self.data.getif(r),self.data.getpol(r), self.data.getcycle(r))459 asaplog.push(out )458 out = " Scan[%d] Beam[%d] IF[%d] Pol[%d] Cycle[%d]" % (scan.getscan(r),scan.getbeam(r),scan.getif(r),scan.getpol(r), scan.getcycle(r)) 459 asaplog.push(out, False) 460 460 self.x = scan._getabcissa(r) 461 461 self.y = scan._getspectrum(r) … … 463 463 self.fit() 464 464 x = self.get_parameters() 465 scan._setspectrum(self.fitter.getresidual(), iRow)465 scan._setspectrum(self.fitter.getresidual(), r) 466 466 print_log() 467 467 return scan -
trunk/python/asaplinefind.py
r522 r880 16 16 print "No lines found!" 17 17 sc2=sc.poly_baseline(fl.get_mask(),7) 18 18 19 19 The algorithm involves a simple threshold criterion. The line is 20 20 considered to be detected if a specified number of consequtive … … 26 26 lowest variances across the spectrum (i.e. histogram equalization is 27 27 used to avoid missing weak lines if strong ones are present). For 28 bad baseline shapes it is rec commended to increase the threshold and28 bad baseline shapes it is recommended to increase the threshold and 29 29 possibly switch the averaging option off (see set_options) to 30 30 detect strong lines only, fit a high order baseline and repeat the line 31 search. 31 search. 32 32 33 33 """ 34 34 35 35 def __init__(self): 36 37 38 39 40 36 """ 37 Create a line finder object. 38 """ 39 self.finder = _asap.linefinder() 40 return 41 41 42 42 def set_options(self,threshold=1.7320508075688772,min_nchan=3, 43 43 avg_limit=8,box_size=0.2): 44 45 46 47 threshold a single channel S/N ratio above which the 48 49 50 51 min_nchan a minimal number of consequtive channels, 44 """ 45 Set the parameters of the algorithm 46 Parameters: 47 threshold a single channel S/N ratio above which the 48 channel is considered to be a detection 49 Default is sqrt(3), which together with 50 min_nchan=3 gives a 3-sigma criterion 51 min_nchan a minimal number of consequtive channels, 52 52 which should satisfy a threshold criterion to 53 54 55 56 57 53 be a detection. Default is 3. 54 avg_limit A number of consequtive channels not greater than 55 this parameter can be averaged to search for 56 broad lines. Default is 8. 57 box_size A running mean box size specified as a fraction 58 58 of the total spectrum length. Default is 1/5 59 Note: For bad baselines threshold should be increased, 60 61 62 undulations instead of real lines. 59 Note: For bad baselines threshold should be increased, 60 and avg_limit decreased (or even switched off completely by 61 setting this parameter to 1) to avoid detecting baseline 62 undulations instead of real lines. 63 63 """ 64 64 self.finder.setoptions(threshold,min_nchan,avg_limit,box_size) 65 66 65 return 66 67 67 def set_scan(self,scan,mask=None,edge=(0,0)): 68 """69 Set the 'data' (scantable) to work with.70 Parameters:71 scan: a scantable72 mask: an optional mask retreived from scantable73 edge: an optional number of channel to drop at74 the edge of spectrum. If only one value is75 specified, the same number will be dropped from76 both sides of the spectrum. Default is to keep77 all channels78 68 """ 79 if not scan: 80 raise RuntimeError, 'Please give a correct scan' 81 82 if isinstance(edge,int): 83 edge=(edge,) 69 Set the 'data' (scantable) to work with. 70 Parameters: 71 scan: a scantable 72 mask: an optional mask retreived from scantable 73 edge: an optional number of channel to drop at 74 the edge of spectrum. If only one value is 75 specified, the same number will be dropped from 76 both sides of the spectrum. Default is to keep 77 all channels 78 """ 79 if not scan: 80 raise RuntimeError, 'Please give a correct scan' 81 if not scan._check_ifs(): 82 raise RuntimeError, 'IFs with different numbers of channels are not yet supported' 83 84 if isinstance(edge,int): 85 edge=(edge,) 84 86 85 87 from asap import _is_sequence_or_number as _is_valid … … 89 91 a pair of integers specified as a tuple" 90 92 91 92 93 if len(edge)>2: 94 raise RuntimeError, "The edge parameter should have two \ 93 95 or less elements" 94 95 96 self.finder.setscan(scan,ones(scan.nchan()),tuple(edge))97 else: 98 99 return 96 if mask is None: 97 from numarray import ones 98 self.finder.setscan(scan,ones(scan.nchan(-1)),tuple(edge)) 99 else: 100 self.finder.setscan(scan,mask,tuple(edge)) 101 return 100 102 def find_lines(self,nRow=0): 101 102 103 104 105 106 103 """ 104 Search for spectral lines in the scan assigned in set_scan. 105 Current Beam/IF/Pol is used, Row is specified by parameter 106 A number of lines found will be returned 107 """ 108 return self.finder.findlines(nRow) 107 109 def get_mask(self,invert=False): 108 109 110 """ 111 Get the mask to mask out all lines that have been found (default) 110 112 111 112 113 Parameters: 114 invert if True, only channels belong to lines will be unmasked 113 115 114 115 116 116 Note: all channels originally masked by the input mask or 117 dropped out by the edge parameter will still be excluded 118 regardless on the invert option 117 119 """ 118 120 return self.finder.getmask(invert) 119 121 def get_ranges(self,defunits=True): 120 121 122 lines found. 122 """ 123 Get ranges (start and end channels or velocities) for all spectral 124 lines found. 123 125 124 125 126 127 if False, the range will be expressed in channels 128 129 130 131 132 126 Parameters: 127 defunits if True (default), the range will use the same units 128 as set for the scan (e.g. LSR velocity) 129 if False, the range will be expressed in channels 130 """ 131 if (defunits): 132 return self.finder.getlineranges() 133 else: 134 return self.finder.getlinerangesinchannels() -
trunk/python/scantable.py
r876 r880 1020 1020 else: return s 1021 1021 1022 # def auto_poly_baseline(self, mask=None, edge=(0,0), order=0, 1023 # threshold=3, insitu=None): 1024 # """ 1025 # Return a scan which has been baselined (all rows) by a polynomial. 1026 # Spectral lines are detected first using linefinder and masked out 1027 # to avoid them affecting the baseline solution. 1028 # 1029 # Parameters: 1030 # mask: an optional mask retreived from scantable 1031 # edge: an optional number of channel to drop at 1032 # the edge of spectrum. If only one value is 1033 # specified, the same number will be dropped from 1034 # both sides of the spectrum. Default is to keep 1035 # all channels 1036 # order: the order of the polynomial (default is 0) 1037 # threshold: the threshold used by line finder. It is better to 1038 # keep it large as only strong lines affect the 1039 # baseline solution. 1040 # insitu: if False a new scantable is returned. 1041 # Otherwise, the scaling is done in-situ 1042 # The default is taken from .asaprc (False) 1043 # 1044 # Example: 1045 # scan2=scan.auto_poly_baseline(order=7) 1046 # """ 1047 # if insitu is None: insitu = rcParams['insitu'] 1048 # varlist = vars() 1049 # from asap.asapfitter import fitter 1050 # from asap.asaplinefind import linefinder 1051 # from asap import _is_sequence_or_number as _is_valid 1052 # 1053 # if not _is_valid(edge, int): 1054 # raise RuntimeError, "Parameter 'edge' has to be an integer or a \ 1055 # pair of integers specified as a tuple" 1056 # 1057 # # setup fitter 1058 # f = fitter() 1059 # f.set_function(poly=order) 1060 # 1061 # # setup line finder 1062 # fl=linefinder() 1063 # fl.set_options(threshold=threshold) 1064 # 1065 # if not insitu: 1066 # workscan=self.copy() 1067 # else: 1068 # workscan=self 1069 # 1070 # rows=range(workscan.nrow()) 1071 # from asap import asaplog 1072 # for i in rows: 1073 # asaplog.push("Processing:") 1074 # asaplog.push(msg) 1075 # fl.set_scan(workscan,mask,edge) 1076 # fl.find_lines(i) 1077 # f.set_scan(workscan, fl.get_mask()) 1078 # f.x = workscan._getabcissa(i) 1079 # f.y = workscan._getspectrum(i) 1080 # f.data = None 1081 # f.fit() 1082 # x = f.get_parameters() 1083 # workscan._setspectrum(f.fitter.getresidual(), i) 1084 # workscan._add_history("poly_baseline", varlist) 1085 # if insitu: 1086 # self._assign(workscan) 1087 # else: 1088 # return workscan 1089 # 1022 def auto_poly_baseline(self, mask=None, edge=(0,0), order=0, 1023 threshold=3, insitu=None): 1024 """ 1025 Return a scan which has been baselined (all rows) by a polynomial. 1026 Spectral lines are detected first using linefinder and masked out 1027 to avoid them affecting the baseline solution. 1028 1029 Parameters: 1030 mask: an optional mask retreived from scantable 1031 edge: an optional number of channel to drop at 1032 the edge of spectrum. If only one value is 1033 specified, the same number will be dropped from 1034 both sides of the spectrum. Default is to keep 1035 all channels 1036 order: the order of the polynomial (default is 0) 1037 threshold: the threshold used by line finder. It is better to 1038 keep it large as only strong lines affect the 1039 baseline solution. 1040 insitu: if False a new scantable is returned. 1041 Otherwise, the scaling is done in-situ 1042 The default is taken from .asaprc (False) 1043 1044 Example: 1045 scan2=scan.auto_poly_baseline(order=7) 1046 """ 1047 if insitu is None: insitu = rcParams['insitu'] 1048 varlist = vars() 1049 from asap.asapfitter import fitter 1050 from asap.asaplinefind import linefinder 1051 from asap import _is_sequence_or_number as _is_valid 1052 1053 if not _is_valid(edge, int): 1054 raise RuntimeError, "Parameter 'edge' has to be an integer or a \ 1055 pair of integers specified as a tuple" 1056 1057 # setup fitter 1058 f = fitter() 1059 f.set_function(poly=order) 1060 1061 # setup line finder 1062 fl=linefinder() 1063 fl.set_options(threshold=threshold) 1064 1065 if not insitu: 1066 workscan=self.copy() 1067 else: 1068 workscan=self 1069 1070 rows=range(workscan.nrow()) 1071 from asap import asaplog 1072 asaplog.push("Processing:") 1073 for r in rows: 1074 msg = " Scan[%d] Beam[%d] IF[%d] Pol[%d] Cycle[%d]" % (workscan.getscan(r),workscan.getbeam(r),workscan.getif(r),workscan.getpol(r), workscan.getcycle(r)) 1075 asaplog.push(msg, False) 1076 fl.set_scan(workscan, mask, edge) 1077 fl.find_lines(r) 1078 f.set_scan(workscan, fl.get_mask()) 1079 f.x = workscan._getabcissa(r) 1080 f.y = workscan._getspectrum(r) 1081 f.data = None 1082 f.fit() 1083 x = f.get_parameters() 1084 workscan._setspectrum(f.fitter.getresidual(), r) 1085 workscan._add_history("poly_baseline", varlist) 1086 if insitu: 1087 self._assign(workscan) 1088 else: 1089 return workscan 1090 1090 1091 # def rotate_linpolphase(self, angle): 1091 1092 # """
Note:
See TracChangeset
for help on using the changeset viewer.