| Index: trunk/tools/wsor/editor_lifecycle/scripts/fitting |
| — | — | @@ -53,6 +53,8 @@ |
| 54 | 54 | parent.add_argument('-minsize', type=int) |
| 55 | 55 | parent.add_argument('-c', '--chop-tail', type=int, default=0, metavar='NUM', |
| 56 | 56 | help='remove %(metavar)s tail observations') |
| | 57 | +parent.add_argument('-xmax', type=int, help='show only (0,xmax) interval') |
| | 58 | +parent.add_argument('-snr', type=float, help='minimum signal-to-noise ratio') |
| 57 | 59 | |
| 58 | 60 | parser = ArgumentParser(description=__doc__) |
| 59 | 61 | subparsers = parser.add_subparsers(help='Parametric models supported') |
| — | — | @@ -80,13 +82,16 @@ |
| 81 | 83 | yy = model(xx) |
| 82 | 84 | kwargs['ecolor'] = 'none' |
| 83 | 85 | kwargs['ls'] = 'none' |
| 84 | | - pp.errorbar(x, y, ye / 2, label=label, **kwargs) |
| | 86 | + kwargs['mfc'] = kwargs.get('color', 'k') |
| | 87 | + kwargs['ms'] = 4 |
| | 88 | + pp.errorbar(x, y, ye, label=label, **kwargs) |
| | 89 | +# pp.plot(x, y, label=label, **kwargs) |
| 85 | 90 | model_label = model.name.split() |
| 86 | 91 | if len(model_label) > 1: |
| 87 | 92 | model_label[1] = model_label[1][:3] + '.' |
| 88 | 93 | model_label = ' '.join(model_label[:2]).capitalize() |
| 89 | 94 | c = kwargs.get('color', 'r') |
| 90 | | - pp.plot(xx, yy, '--', color=c, label='{} fit'.format(model_label), lw=2.5) |
| | 95 | + pp.plot(xx, yy, '--', color=c, lw=2) |
| 91 | 96 | if ns.loglog: |
| 92 | 97 | pp.xscale('log') |
| 93 | 98 | pp.yscale('log') |
| — | — | @@ -105,6 +110,8 @@ |
| 106 | 111 | # r /= np.abs(rm) |
| 107 | 112 | pp.axhline(y=0, c='k') |
| 108 | 113 | kwargs['ls'] = ':' |
| | 114 | + kwargs['lw'] = 2 |
| | 115 | + kwargs['ms'] = 4 |
| 109 | 116 | pp.plot(x, r, label=label, **kwargs) |
| 110 | 117 | pp.title('Fit residuals') |
| 111 | 118 | pp.xlabel('Days since registration') |
| — | — | @@ -143,13 +150,17 @@ |
| 144 | 151 | * both : A = y[0] and C = y[-1] |
| 145 | 152 | ''' |
| 146 | 153 | x, y, ye, n = np.loadtxt(data_file, unpack=True) |
| 147 | | - x = x[:-ns.chop_tail] |
| 148 | | - y = y[:-ns.chop_tail] |
| 149 | | - ye = ye[:-ns.chop_tail] |
| 150 | | - n = n[:-ns.chop_tail] |
| 151 | | - idx = ye > 0 |
| 152 | | - if minsize: |
| 153 | | - idx = idx * (n > minsize) |
| | 154 | + # XXX using the S/N ratio to filter out noisy estimates of the ratio |
| | 155 | + s = ye / np.sqrt(n) |
| | 156 | + snr = y / s |
| | 157 | + idx = snr >= ns.snr |
| | 158 | +# x = x[:-ns.chop_tail] |
| | 159 | +# y = y[:-ns.chop_tail] |
| | 160 | +# ye = ye[:-ns.chop_tail] |
| | 161 | +# n = n[:-ns.chop_tail] |
| | 162 | + idx = idx * (ye > 0) |
| | 163 | +# if minsize: |
| | 164 | +# idx = idx * (n > minsize) |
| 154 | 165 | if not np.any(idx): |
| 155 | 166 | raise EmptyDataError() |
| 156 | 167 | if np.sum(idx) < len(modelclass.__params__): |
| — | — | @@ -157,6 +168,11 @@ |
| 158 | 169 | x = x[idx] |
| 159 | 170 | y = y[idx] |
| 160 | 171 | ye = ye[idx] |
| | 172 | + n = n[idx] |
| | 173 | + # XXX <Sat Oct 29 17:30:38 CEST 2011> getting the standard error of the |
| | 174 | + # mean. This should be really computed in in lifecycle.rates and new data |
| | 175 | + # files be produced using mkrates. |
| | 176 | + ye = ye / np.sqrt(n) |
| 161 | 177 | model = modelclass() |
| 162 | 178 | if constrain in ['head', 'both']: |
| 163 | 179 | model.A = y[np.argmin(np.abs(x))] |
| — | — | @@ -177,11 +193,14 @@ |
| 178 | 194 | |
| 179 | 195 | def main(ns): |
| 180 | 196 | # import matploblib |
| 181 | | - from matplotlib.lines import lineMarkers as markers |
| 182 | | - markers = dict(filter( |
| 183 | | - lambda k : isinstance(k[0],str) and k[1] is not '_draw_nothing', |
| 184 | | - markers.items())).keys() |
| 185 | | - colors = 'bgcmykw' |
| | 197 | +# from matplotlib.lines import lineMarkers as markers |
| | 198 | +# markers = dict(filter( |
| | 199 | +# lambda k : isinstance(k[0],str) and k[1] is not '_draw_nothing', |
| | 200 | +# markers.items())).keys() |
| | 201 | +# markers.remove('.') |
| | 202 | +# markers.remove(',') |
| | 203 | + markers = 'oD^vs*hHpdv<>' |
| | 204 | + colors = 'brgkcmyw' |
| 186 | 205 | |
| 187 | 206 | data = [] |
| 188 | 207 | models = [] |
| — | — | @@ -199,7 +218,8 @@ |
| 200 | 219 | continue |
| 201 | 220 | data.append((x,y,ye)) |
| 202 | 221 | models.append(model) |
| 203 | | - labels.append(os.path.splitext(f)[0]) |
| | 222 | + lab = os.path.splitext(os.path.basename(f))[0] |
| | 223 | + labels.append(lab.replace('_','\\_')) |
| 204 | 224 | |
| 205 | 225 | # plot fits |
| 206 | 226 | pp.figure() |
| — | — | @@ -208,6 +228,8 @@ |
| 209 | 229 | m = markers[i % len(markers)] |
| 210 | 230 | c = colors[i % len(colors)] |
| 211 | 231 | plotfit(model, x, y, ye, label, marker=m, color=c) |
| | 232 | + if ns.xmax is not None: |
| | 233 | + pp.xlim(0, ns.xmax) |
| 212 | 234 | |
| 213 | 235 | # plot residuals |
| 214 | 236 | pp.figure() |
| — | — | @@ -216,6 +238,8 @@ |
| 217 | 239 | m = markers[i % len(markers)] |
| 218 | 240 | c = colors[i % len(colors)] |
| 219 | 241 | plotresid(model, x, y, label, marker=m, color=c) |
| | 242 | + if ns.xmax is not None: |
| | 243 | + pp.xlim(0, ns.xmax) |
| 220 | 244 | |
| 221 | 245 | # save figures |
| 222 | 246 | if ns.output_file is not None: |