1
2
3 """
4 A plotting environment for gas transitions and all that is associated with that.
5
6 Author: R. Lombaert
7
8 """
9
10 import os
11 from scipy import array
12 import operator
13 import subprocess
14 from scipy.interpolate import interp1d
15 import numpy as np
16
17 import cc.path
18 from cc.plotting.objects.PlottingSession import PlottingSession
19 from cc.tools.io import DataIO
20 from cc.modeling.objects import Transition
21 from cc.plotting import Plotting2
22 from cc.tools.readers import LineList
23 from cc.data.instruments import Pacs
24 from cc.modeling.objects import Star
25
26
27
29
30 """
31 Class for plotting gas lines and their information.
32
33 """
34
35 - def __init__(self,star_name,path_gastronoom='Output2014',\
36 inputfilename=None,pacs=None,spire=None,fn_add_star=1):
37
38 """
39 Initializing an instance of PlotGas.
40
41 @param star_name: name of the star from Star.dat, use default only
42 when never using any star model specific things
43 @type star_name: string
44
45 @keyword path_gastronoom: Output modeling folder in MCMax home folder
46
47 (default: 'codeJun2013')
48 @type path_gastronoom: string
49 @keyword inputfilename: name of inputfile that is also copied to the
50 output folder of the plots,
51 if None nothing is copied
52
53 (default: None)
54 @type inputfilename: string
55 @keyword pacs: a Pacs object needed for plotting PACS data. None of no
56 PACS data involved.
57
58 (default: None)
59 @type pacs: Pacs()
60 @keyword spire: a Spire object needed for plotting SPIRE data. None of no
61 SPIRE data involved.
62
63 (default: None)
64 @type spire: Spire()
65 @keyword fn_add_star: Add the star name to the requested plot filename.
66 Only relevant if fn_plt is given in a sub method.
67
68 (default: 1)
69 @type fn_add_star: bool
70
71
72 """
73
74 super(PlotGas, self).__init__(star_name=star_name,\
75 path=path_gastronoom,\
76 code='GASTRoNOoM',\
77 inputfilename=inputfilename,\
78 fn_add_star=fn_add_star)
79
80 cc.path.gout = os.path.join(cc.path.gastronoom,self.path)
81 self.pacs = pacs
82 self.spire = spire
83 self.sphinx_flux_list = []
84
85
86
88
89 '''
90 Make a Star list based on either GASTRoNOoM cooling ids or PACS ids.
91
92 @param models: model_ids for the MCMax db
93 @type models: list(string)
94
95 @return: the parameter sets
96 @rtype: list[Star()]
97
98 '''
99
100 star_grid = Star.makeStars(models=models,\
101 id_type='pacs' in models[0].lower() \
102 and 'PACS' \
103 or 'GASTRoNOoM',\
104 code='GASTRoNOoM',path=self.path)
105 if 'pacs' in models[0].lower():
106 self.pacs.addStarPars(star_grid)
107 [star.addCoolingPars() for star in star_grid]
108 return star_grid
109
110
111
112 - def plotVelocity(self,star_grid=[],models=[],fn_plt='',cfg=''):
113
114 '''
115 Plot velocity versus radius for every model in star_grid.
116
117 @keyword star_grid: List of Star() instances. If default, model ids
118 have to be given.
119
120 (default: [])
121 @type star_grid: list[Star()]
122 @keyword models: The model ids, only required if star_grid is []
123
124 (default: [])
125 @type models: list[string]
126 @keyword fn_plt: A base plot filename. Includes folder. If not, a
127 default is added
128
129 (default: '')
130 @type fn_plt: string
131 @keyword cfg: path to the Plotting2.plotCols config file. If default,
132 the hard-coded default plotting options are used.
133
134 (default: '')
135 @type cfg: string
136
137 '''
138
139 print '***********************************'
140 print '** Plotting Velocity Profiles'
141 if not star_grid and models:
142 star_grid = self.makeStars(models=models)
143 elif (not models and not star_grid) or (models and star_grid):
144 print '** Input is undefined or doubly defined. Aborting.'
145 return
146
147 cfg_dict = Plotting2.readCfg(cfg)
148 if cfg_dict.has_key('filename'):
149 fn_plt = cfg_dict.pop('filename')
150
151 pfns = []
152 for i,star in enumerate(star_grid):
153 if star['LAST_GASTRONOOM_MODEL']:
154 rad = star.getGasRad(unit='rstar')
155 vel = star.getGasVelocity()
156 vel = vel/10.**5
157 avgdrift = star.getAverageDrift()/10.**5
158
159
160 pfn = fn_plt if fn_plt else 'velocity'
161 suff = '{}_{}'.format(star['LAST_GASTRONOOM_MODEL'],i)
162 pfn = self.setFnPlt(pfn,fn_suffix=suff)
163
164 plot_title = '%s %s: Velocity Profile for Model %i (%s)'\
165 %(self.plot_id.replace('_','\_'),\
166 self.star_name_plots,i,\
167 star['LAST_GASTRONOOM_MODEL'].replace('_','\_'))
168 pfns.append(Plotting2.plotCols(x=rad,\
169 y=[vel,avgdrift],cfg=cfg_dict,\
170 filename=pfn,\
171 xaxis='R (R$_*$)',\
172 yaxis=r'$v$ (km s$^{-1}$)',\
173 plot_title=plot_title,\
174 key_location=(0.0,0.0),\
175 keytags=['Gas Velocity',\
176 'Grain-size Weighted Drift'],\
177 xlogscale=1))
178 if pfns and pfns[0][-4:] == '.pdf':
179 pfn = fn_plt if fn_plt else 'velocity_profiles'
180 pfn = self.setFnPlt(pfn) + '.pdf'
181 DataIO.joinPdf(old=sorted(pfns),new=pfn,del_old=0)
182 print '** Plots can be found at:'
183 print pfn
184 print '***********************************'
185 elif pfns:
186 print '** Plots can be found at:'
187 print '\n'.join(pfns)
188 print '***********************************'
189 else:
190 print '** No GASTRoNOoM models were calculated successfully. '+\
191 'No velocity profiles can be plotted.'
192 print '***********************************'
193
194
195
196 - def plotTemp(self,star_grid=[],models=[],fn_plt='',force_plot=0,cfg=''):
197
198 '''
199 Plot temperature profiles of all models.
200
201 @keyword star_grid: List of Star() instances. If default, model ids
202 have to be given.
203
204 (default: [])
205 @type star_grid: list[Star()]
206 @keyword models: The model ids, only required if star_grid is []
207
208 (default: [])
209 @type models: list[string]
210 @keyword fn_plt: A base plot filename. Includes folder. If not, a
211 default is added
212
213 (default: '')
214 @type fn_plt: string
215 @keyword force_plot: force a plotting if more than models are requested
216
217 (default: 0)
218 @type force_plot: bool
219 @keyword cfg: path to the Plotting2.plotCols config file. If default,
220 the hard-coded default plotting options are used.
221
222 (default: '')
223 @type cfg: string
224
225 '''
226
227 print '***********************************'
228 print '** Plotting Gas Temperature Profiles'
229 if not star_grid and models:
230 star_grid = self.makeStars(models=models)
231 elif (not models and not star_grid) or (models and star_grid):
232 print '** Input is undefined or doubly defined. Aborting.'
233 return
234
235 cfg_dict = Plotting2.readCfg(cfg)
236 if cfg_dict.has_key('filename'):
237 fn_plt = cfg_dict.pop('filename')
238 if len(star_grid) < 20 or force_plot:
239 valid_sg = [star
240 for star in star_grid
241 if star['LAST_GASTRONOOM_MODEL']]
242 radii_rstar = [star.getGasRad(unit='rstar') for star in valid_sg]
243 radii = [star.getGasRad(unit='cm') for star in valid_sg]
244 temps = [star.getGasTemperature() for star in valid_sg]
245
246 if temps:
247 keytags = star_grid[0].has_key('LAST_PACS_MODEL') \
248 and ['%s, %s, Mdot = %.2e'\
249 %(star['LAST_GASTRONOOM_MODEL']\
250 .replace('_','\_'),\
251 str(star['LAST_PACS_MODEL'])\
252 .replace('_','\_'),\
253 float(star['MDOT_GAS']))
254 for star in star_grid] \
255 or [star['LAST_GASTRONOOM_MODEL'].replace('_','\_')
256 for star in valid_sg]
257
258
259 pfn = fn_plt if fn_plt else 'temperature_profiles'
260 pfn_rstar = self.setFnPlt(pfn,fn_suffix='rstar')
261 pfn = self.setFnPlt(pfn)
262
263
264 pfn_rstar = Plotting2.plotCols(x=radii_rstar,y=temps,\
265 cfg=cfg,xaxis='R (R$_*$)',\
266 filename=pfn_rstar,yaxis='T (K)',\
267 xlogscale=1,ylogscale=1,keytags=keytags)
268 keys_cm = ['Model %i'%(i+1)
269 for i in xrange(len(star_grid))]
270 pfn = Plotting2.plotCols(x=radii,y=temps,cfg=cfg_dict,\
271 filename=pfn,xaxis='$r$ (cm)',\
272 yaxis='$T_\mathrm{g}$ (K)',\
273 figsize=(12.5,8),fontsize_ticklabels=26,\
274 key_location=(0.05,0.05),xlogscale=1,ylogscale=1,\
275 keytags=keys_cm,fontsize_axis=26,fontsize_key=26)
276 print '** Plots can be found at:'
277 print pfn
278 print pfn_rstar
279 print '***********************************'
280 else:
281 print '** No GASTRoNOoM models were calculated successfully.'+\
282 'No temperature profiles can be plotted.'
283 print '***********************************'
284
285
286
287 - def plotTransitions(self,star_grid,cfg='',no_data=0,vg_factor=3,\
288 telescope_label=1,sort_freq=0,sort_molec=0,\
289 no_models=0,limited_axis_labels=0,date_tag=1,\
290 n_max_models=10,fn_plt='',fn_suffix='',fit_vlsr=1,\
291 plot_intrinsic=0,plot_unresolved=0,cont_subtract=1):
292
293 """
294 Plotting beam convolved line profiles in Tmb for both model and data if
295 available.
296
297 @param star_grid: list of Star objects for which the plotting is done.
298 @type star_grid: list[Star]
299 @keyword cfg: path to the Plotting2.plotCols config file. If default,
300 the hard-coded default plotting options are used.
301
302 (default: '')
303 @type cfg: string
304 @keyword no_data: Don't include the data
305 @type no_data: bool
306 @keyword vg_factor: The factor with which the terminal velocity is
307 multiplied. This determines the xrange of the plots
308 @type vg_factor: float
309 @keyword telescope_label: Include a label showing the telescope name
310
311 (default: 1)
312 @type telescope_label: bool
313 @keyword sort_freq: Sort the lines by frequency rather than wavelength.
314
315 (default: 0)
316 @type sort_freq: bool
317 @keyword sort_molec: Sort the lines by molecule. Can be combined with
318 sort_freq
319
320 (default: 0)
321 @type sort_molec: bool
322 @keyword no_models: Only show data for the resolved lines.
323
324 (default: 0)
325 @type no_models: bool
326 @keyword limited_axis_labels: Remove axis labels not at the left or at
327 the bottom of the tiled plot
328
329 (default: 0)
330 @type limited_axis_labels: bool
331 @keyword date_tag: Add a tag to a plot indicating the date of
332 observation. Only available for non-intrinsic obs.
333
334 (default: 1)
335 @type date_tag: bool
336 @keyword n_max_models: Maximum number of models per tile
337
338 (default: 10)
339 @type n_max_models: bool
340 @keyword fn_plt: A base plot filename. Includes folder. If not, a
341 default is added
342
343 (default: '')
344 @type fn_plt: string
345 @keyword fn_suffix: A suffix that is appended to the filename. For
346 instance, when running the plot command for a
347 best fit subgrid of Star() models as to not
348 overwrite the plot of the full grid.
349
350 (default: '')
351 @type fn_suffix: string
352 @keyword fit_vlsr: Show the models after shifting based on
353 the best_vlsr value, instead of the sphinx output
354
355 (default: 1)
356 @type fit_vlsr: bool
357 @keyword plot_intrinsic: Plot the intrinsic profiles instead of the
358 beam-convolved profiles (for instance when
359 comparing GASTRoNOoM models to LIME models)
360
361 (default: 0)
362 @type plot_intrinsic: bool
363 @keyword plot_unresolved: Plot intrinsic line profiles as well, for
364 unresolved data. The data themselves are not
365 added. By default, this is off as the line
366 profiles do not give you a lot information
367 before convolution with the wavelength
368 resolution.
369
370 (default: 0)
371 @type plot_unresolved: bool
372
373 @keyword cont_subtract: Subtract the continuum from the sphinx line
374 profile
375
376 (default: 1)
377 @type cont_subtract: bool
378
379 """
380
381 print '***********************************'
382 print '** Creating Transition plots.'
383
384 cfg_dict = Plotting2.readCfg(cfg)
385 if cfg_dict.has_key('dimensions'):
386 x_dim = int(cfg_dict['dimensions'][0])
387 y_dim = int(cfg_dict['dimensions'][1])
388 else:
389 x_dim, y_dim = 4,3
390 if cfg_dict.has_key('no_data'):
391 no_data = bool(cfg_dict['no_data'])
392 if cfg_dict.has_key('vg_factor'):
393 vg_factor = float(cfg_dict['vg_factor'])
394 if cfg_dict.has_key('telescope_label'):
395 telescope_label = int(cfg_dict['telescope_label'])
396 if cfg_dict.has_key('sort_freq'):
397 sort_freq = int(cfg_dict['sort_freq'])
398 if cfg_dict.has_key('sort_molec'):
399 sort_molec = int(cfg_dict['sort_molec'])
400 if cfg_dict.has_key('no_models'):
401 no_models = int(cfg_dict['no_models'])
402 if cfg_dict.has_key('limited_axis_labels'):
403 limited_axis_labels = cfg_dict['limited_axis_labels']
404 if cfg_dict.has_key('date_tag'):
405 date_tag = int(cfg_dict['date_tag'])
406 if cfg_dict.has_key('n_max_models'):
407 n_max_models = int(cfg_dict['n_max_models'])
408 if cfg_dict.has_key('plot_intrinsic'):
409 plot_intrinsic = int(cfg_dict['plot_intrinsic'])
410 if cfg_dict.has_key('plot_unresolved'):
411 plot_unresolved = int(cfg_dict['plot_unresolved'])
412 if cfg_dict.has_key('fit_vlsr'):
413 fit_vlsr = int(cfg_dict['fit_vlsr'])
414 if cfg_dict.has_key('cont_subtract'):
415 cont_subtract = int(cfg_dict['cont_subtract'])
416 if fn_plt:
417 if not cfg_dict.has_key('filename'): cfg_dict['filename'] = fn_plt
418 if fn_suffix:
419 filename = cfg_dict.get('filename',None)
420 if not filename is None: filename = '_'.join(filename,fn_suffix)
421 cfg_dict['filename'] = filename
422 if cfg_dict.has_key('keytags'):
423 keytags = cfg_dict['keytags']
424 unreso_keytags = keytags
425 else:
426
427 keytags = [',\\ '.join(set(\
428 [trans.getModelId() != '' \
429 and str(trans.getModelId())\
430 .replace('_','\_')
431 or str(None)
432 for trans in star['GAS_LINES']]))
433 for star in star_grid]
434 unreso_keytags = list(keytags)
435
436
437
438 trans_list = Transition.extractTransFromStars(star_grid,sort_freq,\
439 sort_molec,\
440 dtype='resolved',\
441 reset_data=0)
442
443
444
445 if plot_unresolved:
446 unreso_list = Transition.extractTransFromStars(star_grid,\
447 sort_freq,\
448 sort_molec,\
449 dtype='unresolved')
450 else:
451 unreso_list = []
452
453 def createTilePlots(trans_list,x_dim,y_dim,no_data,intrinsic,\
454 vg_factor,keytags,telescope_label,no_models,cfg,\
455 star_grid,limited_axis_labels,date_tag,indexi,\
456 indexf,fit_vlsr,cont_subtract):
457
458 '''
459 Create a tiled plot for a transition list.
460
461 A list of transitions is exhausted for every star in star_grid,
462 as long as tiles in a single plot are still available.
463
464 @param trans_list: The transition list. Transitions will be removed
465 from this list as tiles are created.
466 @type trans_list: list[Transition()]
467 @param x_dim: The number of tiles in the horizontal direction
468 @type x_dim: int
469 @param y_dim: The number of tiles in the vertical direction
470 @type y_dim: int
471 @param no_data: Include data or not? Will call a function that
472 gathers the data
473 @type no_data: bool
474 @param intrinsic: Intrinsic line profiles, or convolved with beam
475 profile? Set to True for the former, False for
476 the latter
477 @type intrinsic: bool
478 @param star_grid: The grid of models for which the chosen
479 transitions in trans_list are plotted. Can be a
480 subgrid of the original star_grid. (determined
481 before calling this method)
482 @type star_grid: list[Star()]
483 @param keytags: list of keys for the different models in star_grid
484 @type keytags: list[string]
485 @param cfg: The config filename passed to the plotTiles method
486 @type cfg: string/dict
487 @param vg_factor: The factor with which the terminal velocity is
488 multiplied. This determines the xrange of the
489 plots
490 @type vg_factor: float
491 @param telescope_label: Include a label showing the telescope name
492 @type telescope_label: bool
493 @param no_models: Only show data for the resolved lines.
494 @type no_models: bool
495 @param limited_axis_labels: Remove axis labels not at the left or
496 at the bottom of the tiled plot
497
498 (default: 0)
499 @type limited_axis_labels: bool
500 @param date_tag: Add a tag to a plot indicating the date of
501 observation. Only available for non-intrinsic obs.
502 @type date_tag: bool
503 @param indexi: The start index of the models in the star_grid
504 @type indexi: int
505 @param indexf: The end index of the models in the star_grid
506 @type indexf: int
507 @param fit_vlsr: Show the models after shifting based on
508 the best_vlsr value, instead of the sphinx output
509 @type fit_vlsr: bool
510 @param cont_subtract: Subtract the continuum value outside the line
511 from the whole line profile.
512 @type cont_subtract: bool
513
514 @return: The data list with dictionaries for every tile is returned
515 @rtype: list[dict]
516
517 '''
518
519 if cfg.has_key('filename'):
520 fn_plt = cfg.pop('filename')
521 mfn = 'models{}to{}'.format(indexi,indexf)
522 fn_suffix = '' if no_models else mfn
523 else:
524 mfn = '' if no_models else 'models{}to{}'.format(indexi,indexf)
525 fn_suffix = '{}{}'.format(no_data and 'nodata_' or '',mfn)
526 fn_plt = '{}lps'.format(intrinsic and 'intrinsic_' or '')
527
528 missing_trans = 0
529 n_subplots = (x_dim*y_dim) - (keytags and 1 or 0)
530 plot_filenames = []
531 i = 0
532 vexp = max([s['VEL_INFINITY_GAS'] for s in star_grid])
533 while trans_list:
534 i += 1
535 data = []
536
537 ndata = 0
538 for j in xrange(n_subplots):
539 current_trans = trans_list.pop(0)
540 current_sub = [star.getTransition(current_trans)
541 for star in star_grid]
542 if None in current_sub:
543 missing_trans += 1
544
545 if not no_data:
546 current_trans.readData()
547 vlsr = current_trans.getVlsr()
548 noise = current_trans.getNoise()
549 else:
550 vlsr = 0.0
551 noise = None
552 for trans in current_sub:
553 if not trans is None:
554 trans.readSphinx()
555
556
557
558
559
560
561
562
563 trans.setData(current_trans)
564 ddict = dict()
565 ddict['x'] = []
566 ddict['y'] = []
567 if not no_models:
568 for trans in current_sub:
569 if trans is None or trans.sphinx is None:
570 ddict['x'].append([])
571 ddict['y'].append([])
572 continue
573
574 if intrinsic:
575 mvel = trans.sphinx.getVelocityIntrinsic()
576 mlp = trans.sphinx.getLPIntrinsic(cont_subtract\
577 =cont_subtract)
578 mlp = mlp*10**(23)
579 else:
580
581
582
583
584 bvlsr = fit_vlsr and trans.getBestVlsr() or vlsr
585 mvel = trans.sphinx.getVelocity() + bvlsr
586 mlp = trans.sphinx.getLPTmb(cont_subtract=\
587 cont_subtract)
588 ddict['x'].append(mvel)
589 ddict['y'].append(mlp)
590
591
592
593
594 if current_trans.lpdata and not no_data:
595 ddict['histoplot'] = []
596 n_models = len(ddict['x'])
597 for ilp,lp in enumerate(current_trans.lpdata):
598 ddict['x'].append(lp.getVelocity())
599 ddict['y'].append(lp.getFlux())
600 ddict['histoplot'].append(n_models+ilp)
601 if ilp == ndata:
602 ndata += 1
603 ddict['labels'] = \
604 [('%s'%(current_trans.molecule.molecule_plot),0.05,0.87),\
605 ('%s'%(current_trans.makeLabel()),0.05,0.76)]
606 if telescope_label:
607 if True in [trans.sphinx.nans_present
608 for trans in current_sub
609 if (not trans is None \
610 and not trans.sphinx is None)]:
611 telescope_string = '%s*'\
612 %current_trans.telescope.replace('-H2O','')\
613 .replace('-CORRB','')
614 else:
615 telescope_string = '%s'\
616 %current_trans.telescope.replace('-H2O','')\
617 .replace('-CORRB','')
618 ddict['labels'].append((telescope_string,0.73,0.85))
619 if current_trans.lpdata and not no_data and date_tag:
620 ddict['labels'].append(\
621 ('; '.join([lp.getDateObs() \
622 for lp in current_trans.lpdata]),\
623 0.05,0.01))
624
625
626
627 ddict['xmax'] = vlsr + vg_factor * vexp
628 ddict['xmin'] = vlsr - vg_factor * vexp
629 if [yi for yi in ddict['y'] if list(yi)]:
630 ddict['ymax'] = max([max(array(yi)[(array(xi) <= ddict['xmax'])* \
631 (array(xi) >= ddict['xmin'])])
632 for xi,yi in zip(ddict['x'],\
633 ddict['y'])
634 if list(yi)])*1.3
635 ddict['ymin'] = min([min(array(yi)[(array(xi) <= ddict['xmax'])* \
636 (array(xi) >= ddict['xmin'])])
637 for xi,yi in zip(ddict['x'],\
638 ddict['y'])
639 if list(yi)])
640
641 if abs(ddict['ymin']) > 2*ddict['ymax']:
642 ddict['ymax'] = 3*ddict['ymax']
643 if limited_axis_labels:
644 if j%x_dim == 0:
645 ddict['yaxis'] = intrinsic \
646 and r'$F_\nu$ (Jy)' \
647 or '$T_\mathrm{mb}$ (K)'
648 else:
649 ddict['yaxis'] = ''
650 if j in xrange(n_subplots-x_dim,n_subplots):
651 ddict['xaxis'] = r'$v$ (km s$^{-1}$)'
652 else:
653 ddict['xaxis'] = ''
654 data.append(ddict)
655 if not trans_list:
656 break
657
658
659 if fn_suffix:
660 suff = fn_suffix+'_trl{}'.format(i)
661 else:
662 suff = 'trl{}'.format(i)
663 pfn = self.setFnPlt(fn_plt,fn_suffix=suff)
664
665
666 if no_models: these_tags = ['Data '+self.star_name_plots]*ndata
667 else: these_tags = keytags+['Data '+self.star_name_plots]*ndata
668 plot_filenames.append(Plotting2.plotTiles(extension='pdf',\
669 data=data,keytags=these_tags,filename=pfn,\
670 xaxis=r'$v$ (km s$^{-1}$)',fontsize_axis=16,cfg=cfg,\
671 yaxis=intrinsic \
672 and r'$F_\nu$ (Jy)' \
673 or '$T_\mathrm{mb}$ (K)',\
674 fontsize_ticklabels=16,dimensions=(x_dim,y_dim),\
675 fontsize_label=20,linewidth=2))
676 if missing_trans:
677 print 'WARNING! %i requested transitions were '%missing_trans+\
678 'not found for a Star(). Within one CC session, this '+\
679 'should not be the case!'
680 print '** %sine profile plots %scan be found at:'\
681 %(intrinsic and 'Intrinsic l' or 'L',\
682 no_data and 'without data ' or '')
683 if plot_filenames and plot_filenames[0][-4:] == '.pdf':
684 pfn = self.setFnPlt(fn_plt,fn_suffix=fn_suffix)+'.pdf'
685 DataIO.joinPdf(old=plot_filenames,new=pfn)
686 print pfn
687 else:
688 print '\n'.join(plot_filenames)
689 print '***********************************'
690
691
692 if trans_list:
693 j = 0
694 while j < len(star_grid):
695 if plot_intrinsic == 1: no_data = 1
696 i = 0
697 subgrid = []
698 subkeys = []
699 while i < n_max_models and i+j < len(star_grid):
700 subgrid.append(star_grid[i+j])
701 if keytags: subkeys.append(keytags[i+j])
702 i += 1
703
704
705 createTilePlots(trans_list=list(trans_list),\
706 vg_factor=vg_factor,\
707 no_data=no_data,cfg=cfg_dict,star_grid=subgrid,\
708 x_dim=x_dim,y_dim=y_dim,keytags=subkeys,\
709 intrinsic=plot_intrinsic,no_models=no_models,\
710 telescope_label=telescope_label,\
711 limited_axis_labels=limited_axis_labels,\
712 date_tag=date_tag,indexi=j,indexf=j+i-1,\
713 fit_vlsr=fit_vlsr,\
714 cont_subtract=cont_subtract)
715 j += i
716 if unreso_list:
717 j = 0
718 while j < len(star_grid):
719 i = 0
720 subgrid = []
721 subkeys = []
722 while i < n_max_models and i+j < len(star_grid):
723 subgrid.append(star_grid[i+j])
724 if keytags: subkeys.append(unreso_keytags[i+j])
725 i += 1
726 createTilePlots(trans_list=list(unreso_list),\
727 vg_factor=vg_factor,\
728 no_data=1,cfg=cfg_dict,star_grid=subgrid,\
729 intrinsic=1,keytags=subkeys,\
730 x_dim=x_dim,y_dim=y_dim,date_tag=date_tag,\
731 telescope_label=telescope_label,no_models=0,\
732 limited_axis_labels=limited_axis_labels,\
733 indexi=j,indexf=j+i-1,fit_vlsr=fit_vlsr,\
734 cont_subtract=cont_subtract)
735 j += i
736
737
738
741
742 '''
743 Create a list of line labels for all molecules and transitions
744 in the molecular linelists requested.
745
746 This is used for spectroscopic databases only! Such as JPL, CDMS, LAMDA
747
748 @param star: The parameter set
749 @type star: Star()
750 @param xmin: minimum wavelength
751 @type xmin: float
752 @param xmax: maximum wavelength
753 @type xmax: float
754
755 @keyword xunit: The unit of the xmax/xmin
756
757 (default: micron)
758 @type xunit: string
759 @keyword fn_trans_marker: A file that includes TRANSITION definitions.
760 These transitions will be marked up in the
761 plot. For instance, when indicating a subset
762 of transitions for one reason or another.
763 The line type can be set for this specific
764 subset, differently from other lines and
765 regardless of the molecule. In combination
766 with a doubly defined line label (through the
767 star_grid['GAS_LINES']), lines can be marked
768 up.
769
770 (default: '')
771 @type fn_trans_marker: string
772 @keyword instrument: The instrument object for which the line labels
773 are created. Used to retrieve the v_lsr. Either
774 'PACS' or 'SPIRE'
775
776 (default: 'PACS')
777 @type instrument: str
778
779 @return: The labels with x location and a molecule index.
780 @rtype: list[string, float, index]
781
782 '''
783
784 cats = star['LL_CAT']
785 min_strengths = star['LL_MIN_STRENGTH']
786 max_excs = star['LL_MAX_EXC']
787 molecs = star['LL_GAS_LIST']
788
789 linelist = []
790 for m,cat,min_str,max_exc in zip(molecs,cats,min_strengths,max_excs):
791
792 if not 'p1H' in m.molecule:
793 fn = os.path.join(cc.path.ll,'{}_{}.dat'.format(m.molecule,cat))
794 ll = LineList.LineList(fn=fn,x_min=xmin,\
795 unit=xunit,x_max=xmax,\
796 min_strength=min_strength,\
797 max_exc=max_exc)
798 linelist.append(ll.makeTransitions(m))
799 lls = self.createLineLabels(linelist=linelist,\
800 fn_trans_marker=fn_trans_marker,
801 instrument=instrument)
802 return lls
803
804
805
806 - def plotLineLists(self,star_grid,include_sphinx=1,cfg='',fn_plt='',\
807 fn_trans_marker='',instrument='PACS'):
808
809 '''
810 Plot linelists along with the indicated data.
811
812 @param star_grid: The Parameter sets
813 @type star_grid: list[Star()]
814
815 @keyword include_sphinx: Include convolved Sphinx models in the plots
816 for the star_grid
817
818 (default: 1)
819 @type include_sphinx: bool
820 @keyword cfg: path to the Plotting2.plotCols config file. If default, the
821 hard-coded default plotting options are used.
822
823 (default: '')
824 @type cfg: string
825 @keyword fn_plt: A base plot filename. Includes folder. If not, a
826 default is added
827
828 (default: '')
829 @type fn_plt: string
830 @keyword fn_trans_marker: A file that includes TRANSITION definitions.
831 These transitions will be marked up in the
832 plot. For instance, when indicating a subset
833 of transitions for one reason or another.
834 The line type can be set for this specific
835 subset, differently from other lines and
836 regardless of the molecule. In combination
837 with a doubly defined line label (through the
838 star_grid['GAS_LINES']), lines can be marked
839 up.
840
841 (default: '')
842 @type fn_trans_marker: string
843 @keyword instrument: The unresolved-data instrument for which the line
844 lists are to be plotted. Either 'PACS' or 'SPIRE'
845
846 (default: 'PACS')
847 @type instrument: str
848
849 '''
850
851 print '***********************************'
852 print '** Starting to plot line identifications for %s from databases.'\
853 %self.star_name
854 cfg_dict = Plotting2.readCfg(cfg)
855 if cfg_dict.has_key('filename'):
856 fn_plt = cfg_dict.pop('filename')
857 if cfg_dict.has_key('instrument'):
858 instrument = cfg_dict['instrument']
859 instrument = instrument.upper()
860 if instrument == 'SPIRE':
861 instr = self.spire
862 else:
863 instrument = 'PACS'
864 instr = self.pacs
865
866 if instr is None:
867 print '** No %s_PATH given. Cannot plot line lists '%instrument + \
868 'without data information. Aborting...'
869 return
870 if cfg_dict.has_key('fn_trans_marker'):
871 fn_trans_marker = cfg_dict['fn_trans_marker']
872 if cfg_dict.has_key('include_sphinx'):
873 include_sphinx = bool(cfg_dict['include_sphinx'])
874 xmins = [min(wave_list) for wave_list in instr.data_wave_list]
875 xmaxs = [max(wave_list) for wave_list in instr.data_wave_list]
876 lls = self.createLineLabelsFromLineLists(star=star_grid[0],\
877 xmin=min(xmins),\
878 xmax=max(xmaxs),\
879 fn_trans_marker=\
880 fn_trans_marker,\
881 instrument=instrument)
882 plot_filenames = []
883 if include_sphinx:
884 if set([s['MOLECULE'] and 1 or 0 for s in star_grid]) \
885 == set([0]) \
886 or set([s['LAST_GASTRONOOM_MODEL'] for s in star_grid]) \
887 == set(['']):
888 include_sphinx = 0
889 else:
890 if instrument == 'PACS':
891 self.setSphinxPacs(star_grid)
892 for i_file,(wave,flux,filename,xmin,xmax) in enumerate(\
893 zip(instr.data_wave_list,instr.data_flux_list,\
894 instr.data_filenames,xmins,xmaxs)):
895 if include_sphinx and instrument == 'PACS':
896 sphinx_flux = [sphinx
897 for sphinx in self.sphinx_flux_list[i_file]
898 if list(sphinx)]
899 elif include_sphinx and instrument == 'SPIRE':
900 sphinx_flux = [instr.getSphinxConvolution(star,filename)[1]
901 for star in star_grid]
902 else:
903 sphinx_flux = []
904
905
906 pfn = fn_plt if fn_plt else 'line_id'
907 suff = os.path.split(filename)[1].replace('.dat','')
908 pfn = self.setFnPlt(pfn,fn_suffix=suff,fn_subfolder='LineLists')
909
910 keytags = ['%s %s'%(instrument,filename.replace('_','\_'))] + \
911 ['Model %i: %s'\
912 %(i+1,instrument=='PACS' \
913 and str(star['LAST_PACS_MODEL']).replace('_','\_')\
914 or star['LAST_GASTRONOOM_MODEL'].replace('_','\_'))
915 for i,star in enumerate(star_grid)
916 if star['LAST_GASTRONOOM_MODEL'] and include_sphinx]
917 plot_filenames.append(Plotting2.plotCols(\
918 x=[wave]*(len(sphinx_flux)+1),y=[flux]+sphinx_flux,\
919 cfg=cfg_dict,filename=pfn,keytags=keytags,\
920 plot_title=self.star_name_plots,histoplot=[0],\
921 number_subplots=3,line_labels=lls,\
922 line_label_color=1,line_label_lines=1,\
923 line_label_spectrum=1))
924
925 pfn = fn_plt if fn_plt else 'line_id'
926 suff = instrument.lower()
927 pfn = self.setFnPlt(pfn,fn_suffix=suff,fn_subfolder='LineLists')
928 pfn += '.pdf'
929 DataIO.joinPdf(old=sorted(plot_filenames),new=pfn)
930 print '** Plots can be found at:'
931 print pfn
932 print '***********************************'
933
934
935
938
939 '''
940 Plot abundance profiles for all molecules in every model.
941
942 @keyword star_grid: List of Star() instances. If default, model ids
943 have to be given.
944
945 (default: [])
946 @type star_grid: list[Star()]
947 @keyword models: The model ids, only required if star_grid is []
948
949 (default: [])
950 @type models: list[string]
951 @keyword cfg: path to the Plotting2.plotCols config file. If default,
952 the hard-coded default plotting options are used.
953
954 (default: '')
955 @type cfg: string
956 @keyword fn_plt: A base plot filename. Includes folder. If not, a
957 default is added
958
959 (default: '')
960 @type fn_plt: string
961 @keyword per_molecule: Plot one molecule for all models in one figure.
962
963 (default: 0)
964 @type per_molecule: bool
965 @keyword unit: The radial unit. Can be 'cm', 'au', 'm' or 'rstar'
966
967 (default: cm)
968 @type unit: str
969
970 '''
971
972 print '***********************************'
973 print '** Plotting Abundance Profiles'
974 if not star_grid and models:
975 star_grid = self.makeStars(models=models)
976 elif (not models and not star_grid) or (models and star_grid):
977 print '** Input is undefined or doubly defined. Aborting.'
978 return
979 pfns = []
980 cfg_dict = Plotting2.readCfg(cfg)
981 if cfg_dict.has_key('filename'):
982 fn_plt = cfg_dict.pop('filename')
983 if cfg_dict.has_key('per_molecule'):
984 per_molecule = cfg_dict['per_molecule']
985 if cfg_dict.has_key('unit'):
986 unit = cfg_dict['unit'].lower()
987
988
989 extra_pars = dict()
990 extra_pars['ymin'] = 1e-9
991 extra_pars['ymax'] = 1e-3
992 extra_pars['ylogscale'] = 1
993 extra_pars['xlogscale'] = 1
994 extra_pars['figsize'] = (12.5,8.5)
995
996 if unit == 'cm': xaxis = '$r\ \mathrm{(cm)}$'
997 elif unit =='au': xaxis = '$r\ \mathrm{(AU)}$'
998 elif unit == 'm': xaxis = '$r\ \mathrm{(m)}$'
999 else: xaxis = '$r\ \mathrm{(R}_\star\mathrm{)}$'
1000 extra_pars['xaxis'] = xaxis
1001
1002
1003 ddata = dict()
1004 for istar,star in enumerate(star_grid):
1005 if not star['LAST_GASTRONOOM_MODEL']: continue
1006 ddata[istar] = dict()
1007 for molec in star['GAS_LIST']:
1008 mid = molec.getModelId()
1009 if not mid: continue
1010 ddata[istar][molec.molecule] = dict()
1011 rad = star.getGasRad(unit=unit,ftype='1',mstr=molec.molecule,
1012 modelid=mid)
1013 nh2 = star.getGasNumberDensity(ftype='1',mstr=molec.molecule,
1014 modelid=mid)
1015 nmol = star.getGasNumberDensity(molecule=1,ftype='1',\
1016 mstr=molec.molecule,\
1017 modelid=mid)
1018 ddata[istar][molec.molecule]['rad'] = rad
1019 ddata[istar][molec.molecule]['nh2'] = nh2
1020 ddata[istar][molec.molecule]['nmol'] = nmol
1021 if molec.set_keyword_change_abundance:
1022 cff = DataIO.readCols(molec.change_fraction_filename)
1023 rfrac,frac = cff[0],cff[1]
1024 rfrac = rfrac
1025 frac_interpol = interp1d(rfrac,frac)(rad)
1026
1027
1028
1029
1030
1031
1032 else:
1033 frac_interpol = 1
1034
1035
1036 abun = nmol/nh2*frac_interpol/molec.abun_factor
1037 ddata[istar][molec.molecule]['abun'] = abun
1038 ddata[istar][molec.molecule]['key'] = molec.molecule_plot
1039 ddata[istar][molec.molecule]['id'] = mid
1040
1041 if not per_molecule:
1042
1043 radii = [dmol['rad'] for molec,dmol in ddata[istar].items()]
1044 abuns = [dmol['abun'] for molec,dmol in ddata[istar].items()]
1045 keytags = [dmol['key'] for molec,dmol in ddata[istar].items()]
1046 ids = [dmol['id'] for molec,dmol in ddata[istar].items()]
1047
1048
1049 if star.has_key('R_DES_H2O'):
1050 radii.extend([array([star['R_DES_H2O'],star['R_DES_H2O']])])
1051 abuns.extend([[1e-2,1e-9]])
1052 keytags.append('Condensation radius H$_2$O ice')
1053 lt.append('--k')
1054 if star['R_OH1612']:
1055 radii.extend([array([star['R_OH1612'],star['R_OH1612']])])
1056 abuns.extend([[1e-2,1e-9]])
1057 keytags.append('Location OH maser')
1058 lt.append('-k')
1059
1060
1061 yaxis = '$n_\mathrm{molec}/n_{\mathrm{H}_2}$'
1062
1063
1064 pfn = fn_plt if fn_plt else 'abundance_profiles'
1065 suff = '_'.join(list(set(ids)))
1066 pfn = self.setFnPlt(pfn,fn_suffix=suff)
1067
1068 pfns.append(Plotting2.plotCols(x=radii,y=abuns,cfg=cfg_dict,\
1069 filename=pfn,keytags=keytags,\
1070 yaxis=yaxis,**extra_pars))
1071
1072 if per_molecule:
1073
1074 molecs = list(set([molec for istar in ddata.keys()
1075 for molec in ddata[istar].keys()]))
1076 for molec in molecs:
1077
1078 mplot = ddata[0][molec]['key']
1079 radii = [dmol['rad']
1080 for istar in ddata.keys()
1081 for imolec,dmol in ddata[istar].items()
1082 if molec == imolec]
1083 abuns = [dmol['abun']
1084 for istar in ddata.keys()
1085 for imolec,dmol in ddata[istar].items()
1086 if molec == imolec]
1087 keytags = [dmol['id'].replace('_','\_')
1088 for istar in ddata.keys()
1089 for imolec,dmol in ddata[istar].items()
1090 if molec == imolec]
1091
1092
1093 yaxis = '$n_\mathrm{%s}/n_{\mathrm{H}_2}$'%mplot.replace('$','')
1094
1095
1096 pfn = fn_plt if fn_plt else 'abundance_profiles'
1097 pfn = self.setFnPlt(pfn,fn_suffix=molec)
1098
1099 pfns.append(Plotting2.plotCols(x=radii,y=abuns,yaxis=yaxis,\
1100 filename=pfn,keytags=keytags,\
1101 cfg=cfg_dict,**extra_pars))
1102
1103 if not per_molecule and pfns and pfns[0][-4:] == '.pdf':
1104 pfn = fn_plt if fn_plt else 'abundance_profiles'
1105 pfn = self.setFnPlt(pfn) + '.pdf'
1106 DataIO.joinPdf(old=pfns,new=pfn)
1107 print '** Plots can be found at:'
1108 print pfn
1109 print '***********************************'
1110 else:
1111 print '** Plots can be found at:'
1112 print '\n'.join(pfns)
1113 print '***********************************'
1114
1115
1116
1117 - def plotLineContributions(self,star_grid,fn_plt='',normalized=1,cfg='',\
1118 do_sort=1,include='intensity'):
1119
1120 '''
1121 Plot the source function as function of impact parameter for every
1122 transition.
1123
1124 @param star_grid: The model parameter sets
1125 @type star_grid: list[Star()]
1126
1127 @keyword fn_plt: A base plot filename. Includes folder. If not, a
1128 default is added
1129
1130 (default: '')
1131 @type fn_plt: string
1132 @keyword cfg: path to the Plotting2.plotCols config file. If default,
1133 the hard-coded default plotting options are used.
1134
1135 (default: '')
1136 @type cfg: string
1137 @keyword normalized: plot the normalized source functions as opposed
1138 to not normalized
1139
1140 (default: 1)
1141 @type normalized: bool
1142 @keyword do_sort: Sort the transition list according to wavelength. If
1143 off, the original order given in the CC input file is
1144 kept
1145
1146 (default: 1)
1147 @type do_sort: bool
1148 @keyword include: Include and additional profile on a second axis. Can
1149 be 'velocity' or 'intensity' (at line center) at the
1150 moment. Any other value will add no second axis.
1151
1152 (default: 'intensity')
1153 @type include: str
1154
1155 '''
1156
1157 print '***********************************'
1158 print '** Plotting Line Contributions'
1159 cfg_dict = Plotting2.readCfg(cfg)
1160 if cfg_dict.has_key('filename'):
1161 fn_plt = cfg_dict.pop('filename')
1162 if cfg_dict.has_key('do_sort'):
1163 do_sort = int(cfg_dict['do_sort'])
1164 if cfg_dict.has_key('normalized'):
1165 normalized = int(cfg_dict['normalized'])
1166 if cfg_dict.has_key('include'):
1167 include = cfg_dict['include']
1168
1169 normalized = int(normalized)
1170 lcf = 'getNormalizedIntensity' if normalized else 'getWeightedIntensity'
1171 for i,star in enumerate(star_grid):
1172 extra_pars = dict()
1173 if do_sort:
1174 transitions = sorted([trans
1175 for trans in star['GAS_LINES']
1176 if trans.getModelId()],\
1177 key=lambda x:x.wavelength)
1178 else:
1179 transitions = star['GAS_LINES']
1180
1181
1182 [t.readSphinx() for t in transitions]
1183 radii = [t.sphinx.getImpact() for t in transitions]
1184 linecontribs = [getattr(t.sphinx,lcf)() for t in transitions]
1185
1186
1187 pfn = fn_plt if fn_plt else 'linecontrib'
1188 subf = 'LCs'
1189 suff = '{}_{}'.format(star['LAST_GASTRONOOM_MODEL'],i)
1190 pfn = self.setFnPlt(pfn,fn_suffix=suff,fn_subfolder=subf)
1191 extra_pars['filename'] = pfn
1192
1193
1194 extra_pars['keytags'] = [r'$I(p,LC)\ pdp$ - $\mathrm{%s}:$ %s'\
1195 %(t.molecule.molecule_plot\
1196 .replace('$',''),\
1197 t.makeLabel())
1198 for t in transitions]
1199 extra_pars['key_location'] = 'best'
1200 extra_pars['ymin'] = normalized and -0.01 or None
1201 extra_pars['ymax'] = normalized and 1.02 or None
1202 extra_pars['xmin'] = 1
1203 extra_pars['xaxis'] = '$p\ \mathrm{(R}_\star\mathrm{)}$'
1204 extra_pars['yaxis'] = '$I(p)\ pdp$'
1205 extra_pars['linewidth'] = 3
1206 extra_pars['xlogscale'] = 1
1207 extra_pars['fontsize_key'] = 16
1208
1209
1210 if include == 'velocity':
1211 rad = star.getGasRad(unit='rstar')
1212 vel = star.getGasVelocity()
1213 vel = vel/10.**5
1214 extra_pars['twiny_x'] = [rad]
1215 extra_pars['twiny_y'] = [vel]
1216 extra_pars['twiny_keytags'] = [r'$v_\mathrm{g}$']
1217 extra_pars['twinyaxis'] = r'$v_\mathrm{g}$ (km s$^{-1}$)'
1218
1219
1220 elif include == 'intensity':
1221 extra_pars['twiny_x'] = radii
1222 extra_pars['twiny_y'] = [t.sphinx.getLineIntensity()
1223 for t in transitions]
1224 extra_pars['twiny_keytags'] = [k.replace(r'\ pdp','')
1225 for k in extra_pars['keytags']]
1226 extra_pars['twiny_logscale'] = 1
1227 extra_pars['twinyaxis'] = r'$I(p)$'
1228
1229 pfn = Plotting2.plotCols(x=radii,y=linecontribs,cfg=cfg_dict,\
1230 **extra_pars)
1231 print '** Plot can be found at:'
1232 print pfn
1233 print '***********************************'
1234
1235
1236
1238
1239 '''
1240 Prepare Sphinx output in Pacs format (including convolution).
1241
1242 @param star_grid: Parameter sets. If empty list, no sphinx models are
1243 set, but an empty list is set for each datafile.
1244 @type star_grid: list[Star()]
1245
1246 @keyword refresh_sphinx_flux: redo the sphinx flux list by pulling from
1247 db, regardless of whether it's already
1248 been done or not.
1249
1250 (default: 0)
1251 @type refresh_sphinx_flux: bool
1252
1253 '''
1254
1255 if not self.sphinx_flux_list or refresh_sphinx_flux:
1256 self.pacs.prepareSphinx(star_grid)
1257
1258
1259 self.sphinx_flux_list = [[self.pacs.getSphinxConvolution(star,fn)[1]
1260 for star in star_grid]
1261 for fn in self.pacs.data_filenames]
1262
1263
1264
1265
1266 - def createLineLabels(self,star_grid=[],linelist=[],molecules=[],\
1267 fn_trans_marker='',\
1268 unit='micron',mark_undetected=0,instrument='PACS'):
1269
1270 '''
1271 Create line labels for all transitions in Star() objects or in
1272 LineList() objects or in a TRANSITION definition file. Priority:
1273 star_grid > linelists. fn_trans_marker is always added in addition.
1274
1275 @keyword star_grid: The Star() models.
1276
1277 (default: [])
1278 @type star_grid: list[Star()]
1279 @keyword linelist: The list of Transition() objects extracted from a
1280 catalog, eg by createLineLabelsFromLineLists.
1281
1282 (default: [])
1283 @type linelist: list[LineList()]
1284 @keyword fn_trans_marker: A file that includes TRANSITION definitions.
1285 These transitions will be marked up in the
1286 plot. For instance, when indicating a subset
1287 of transitions for one reason or another.
1288 The line type can be set for this specific
1289 subset, differently from other lines and
1290 regardless of the molecule. In combination
1291 with a doubly defined line label (through
1292 star_grid/linelists), lines can be marked
1293 up.
1294
1295 (default: '')
1296 @type fn_trans_marker: string
1297 @keyword mark_undetected: Mark the undetected transitions in the same
1298 way extra marked transitions would be marked
1299 by fn_trans_marker.
1300
1301 (default: 0)
1302 @type mark_undetected: bool
1303 @keyword unit: The unit of the location number. Can be 'micron' or
1304 'cm-1'.
1305
1306 (default: 'micron')
1307 @type unit: string
1308 @keyword instrument: The instrument object for which the line labels
1309 are created. Used to retrieve the v_lsr. Either
1310 'PACS' or 'SPIRE'
1311
1312 (default: 'SPIRE')
1313 @type instrument: str
1314
1315 @return: a sorted list(set) of line labels
1316 @rtype: list[string]
1317
1318 '''
1319
1320 instrument = instrument.upper()
1321 if instrument == 'PACS':
1322 vlsr = self.pacs.vlsr
1323 elif instrument == 'SPIRE':
1324 vlsr = self.spire.vlsr
1325
1326 if star_grid:
1327 alltrans = Transition.extractTransFromStars(star_grid,\
1328 dtype=instrument,\
1329 reset_data=0)
1330 elif linelist:
1331 alltrans = linelist
1332 else:
1333 alltrans = []
1334
1335 lls = [('%s %s'%(t.molecule.molecule,t.makeLabel()),\
1336 t.wavelength*10**4*1./(1-vlsr/t.c),\
1337 t.molecule.molecule_index,\
1338 t.vup>0)
1339 for t in alltrans]
1340
1341 used_indices = list(set([ll[-2] for ll in lls]))
1342 if fn_trans_marker:
1343 all_molecs = set([t.molecule for t in alltrans])
1344 def_molecs = dict([(m.molecule,m) for m in all_molecs])
1345 if star_grid: star = star_grid[0]
1346 else: star = None
1347 trl = DataIO.readDict(fn_trans_marker,multi_keys=['TRANSITION'])
1348 n_entry = len(trl['TRANSITION'][0].split())
1349 trl_sorted = DataIO.checkEntryInfo(trl['TRANSITION'],n_entry,\
1350 'TRANSITION')
1351 etrans = [Transition.makeTransition(trans=t,def_molecs=def_molecs,\
1352 star=star)
1353 for t in trl_sorted]
1354 this_index = max(used_indices)+1
1355 used_indices = used_indices + [this_index]
1356 ells = [('%s %s'%(t.molecule.molecule,t.makeLabel()),\
1357 t.wavelength*10**4*1./(1-vlsr/t.c),\
1358 this_index,\
1359 t.vup>0)
1360 for t in etrans]
1361 lls = lls + ells
1362
1363 if mark_undetected:
1364 this_index = max(used_indices)+1
1365 extra_trans = [t
1366 for t in alltrans
1367 if t.getIntIntUnresolved()[0] is None]
1368 ells = [('%s %s'%(t.molecule.molecule,t.makeLabel()),\
1369 t.wavelength*10**4*1./(1-vlsr/t.c),\
1370 this_index,\
1371 t.vup>0)
1372 for t in extra_trans]
1373 lls = lls + ells
1374
1375 if unit == 'cm-1':
1376 lls = [(l,1./w*10**4,i,vib) for l,w,i,vib in lls]
1377 lls = sorted(lls,key=operator.itemgetter(1))
1378 return lls
1379
1380
1381
1382 - def plotPacsLineScans(self,star_grid=[],models=[],exclude_data=0,cfg='',\
1383 cont_subtracted=1,fn_trans_marker='',fn_plt='',\
1384 dimensions=(5,2),mark_undetected=0,\
1385 remove_axis_titles=1,include_band=1):
1386
1387 '''
1388 Plot PACS line scans.
1389
1390 Data can be in- or excluded, as can models.
1391
1392 Both continuum subtracted data as well as the original spectra can be
1393 plotted.
1394
1395 @keyword star_grid: star models for which PACS data will be fetched.
1396
1397 (default: [])
1398 @type star_grid: list[Star()]
1399 @keyword models: list of pacs_ids or gastronoom model ids. If neither
1400 this or star_grid are defined, only data are plotted.
1401 If star_grid is defined, this keyword is ignored.
1402 (default: [])
1403 @type models: list[strings]
1404 @keyword cfg: path to the Plotting2.plotCols config file. If default,
1405 the hard-coded default plotting options are used.
1406
1407 (default: '')
1408 @type cfg: string
1409 @keyword fn_plt: A plot filename for the tiled plot.
1410
1411 (default: '')
1412 @type fn_plt: string
1413 @keyword fn_trans_marker: A file that includes TRANSITION definitions.
1414 These transitions will be marked up in the
1415 plot. For instance, when indicating a subset
1416 of transitions for one reason or another.
1417 The line type can be set for this specific
1418 subset, differently from other lines and
1419 regardless of the molecule. In combination
1420 with a doubly defined line label (through the
1421 star_grid['GAS_LINES']), lines can be marked
1422 up.
1423
1424 (default: '')
1425 @type fn_trans_marker: string
1426 @keyword mark_undetected: Mark the undetected transitions in the same
1427 way extra marked transitions would be marked
1428 by fn_trans_marker.
1429
1430 (default: 0)
1431 @type mark_undetected: bool
1432 @keyword remove_axis_titles: Remove axis titles in between different
1433 tiles in the plot. Only keeps the ones on
1434 the left and the bottom of the full plot.
1435
1436 (default: 1)
1437 @type remove_axis_titles: bool
1438 @keyword include_band: Include a label that names the band.
1439
1440 (default: 1)
1441 @type include_band: bool
1442 @keyword exclude_data: if enabled only the sphinx mdels are plotted.
1443
1444 (default: 0)
1445 @type exclude_data: bool
1446 @keyword cont_subtracted: Plot the continuum subtracted data.
1447
1448 (default: 1)
1449 @type cont_subtracted: bool
1450 @keyword dimensions: The number of tiles in the x and y direction is
1451 given: (x-dim,y-dim)
1452
1453 (default: (5,2))
1454 @type dimensions: tuple(int,int)
1455
1456 '''
1457
1458 print '***********************************'
1459 print '** Plotting line scans.'
1460 if self.pacs is None:
1461 print '** No PATH_PACS given. Cannot plot PACS spectra without '+\
1462 'data information. Aborting...'
1463 return
1464 if not star_grid and models:
1465 star_grid = self.makeStars(models=models)
1466
1467 self.setSphinxPacs(star_grid)
1468 cfg_dict = Plotting2.readCfg(cfg)
1469 if cfg_dict.has_key('filename'):
1470 fn_plt = cfg_dict.pop('filename')
1471 if cfg_dict.has_key('exclude_data'):
1472 exclude_data = bool(cfg_dict['exclude_data'])
1473 if cfg_dict.has_key('fn_trans_marker'):
1474 fn_trans_marker = cfg_dict['fn_trans_marker']
1475 if cfg_dict.has_key('cont_subtracted'):
1476 cont_subtracted = cfg_dict['cont_subtracted']
1477 if cfg_dict.has_key('mark_undetected'):
1478 mark_undetected = cfg_dict['mark_undetected']
1479 if cfg_dict.has_key('remove_axis_titles'):
1480 remove_axis_titles = cfg_dict['remove_axis_titles']
1481 if cfg_dict.has_key('dimensions'):
1482 dimensions = (int(cfg_dict['dimensions'][0]),\
1483 int(cfg_dict['dimensions'][1]))
1484
1485 if not star_grid:
1486 exclude_data = 0
1487
1488 if mark_undetected:
1489 trl = Transition.extractTransFromStars(star_grid,dtype='PACS')
1490 for ifn in range(len(self.pacs.data_filenames)):
1491 self.pacs.intIntMatch(trans_list=trl,ifn=ifn)
1492
1493 lls = self.createLineLabels(star_grid=star_grid,\
1494 fn_trans_marker=fn_trans_marker,\
1495 mark_undetected=mark_undetected,\
1496 instrument='PACS')
1497 tiles = []
1498 print '** Plotting now...'
1499 for idd,(wave,flux,flux_ori,sphinx_flux,filename,ordername) in \
1500 enumerate(sorted(zip(self.pacs.data_wave_list,\
1501 self.pacs.data_flux_list,\
1502 self.pacs.data_original_list,\
1503 self.sphinx_flux_list,\
1504 self.pacs.data_filenames,\
1505 self.pacs.data_ordernames),\
1506 key=lambda x: x[0][0])):
1507 ddict = dict()
1508 ddict['x'] = exclude_data \
1509 and [wave]*(len(sphinx_flux)) \
1510 or [wave]*(len(sphinx_flux)+1)
1511 d_yvals = cont_subtracted and [flux] or [flux_ori]
1512 ddict['y'] = exclude_data and sphinx_flux or d_yvals+sphinx_flux
1513 if include_band:
1514 ddict['labels'] = [(ordername,0.08,0.80)]
1515 ddict['xmin'] = wave[0]
1516 ddict['xmax'] = wave[-1]
1517 ddict['ymin'] = 0.95*min([min(ff) for ff in ddict['y'] if ff.size])
1518 ddict['ymax'] = 1.2*max([max(ff) for ff in ddict['y'] if ff.size])
1519 ddict['histoplot'] = (not exclude_data) and [0] or []
1520 ddict['line_labels'] = lls
1521
1522
1523
1524 if remove_axis_titles:
1525 n_tiles = len(self.pacs.data_filenames)
1526 ddict['xaxis'] = idd in range(n_tiles-dimensions[0],n_tiles)\
1527 and r'$\lambda$ ($\mu$m)' or ''
1528 ddict['yaxis'] = idd%dimensions[0] == 0 \
1529 and r'$F_\nu$ (Jy)' or ''
1530 tiles.append(ddict)
1531
1532
1533 pfn = fn_plt if fn_plt else 'PACS_linescans'
1534 subf = 'PACS_results'
1535 pfn = self.setFnPlt(pfn,fn_subfolder=subf)
1536
1537 pfn = Plotting2.plotTiles(filename=fn_plt,data=tiles,cfg=cfg_dict,\
1538 line_label_color=1,fontsize_label=15,\
1539 line_label_lines=1,dimensions=dimensions)
1540 print '** Your plot can be found at:'
1541 print pfn
1542 print '***********************************'
1543
1544
1545
1546 - def plotPacs(self,star_grid=[],models=[],exclude_data=0,fn_plt='',cfg='',\
1547 fn_trans_marker='',include_band=1,number_subplots=3,\
1548 mark_undetected=0):
1549
1550 '''
1551 Plot PACS data along with Sphinx results, one plot per band.
1552
1553 @keyword star_grid: star models for which PACS data will be fetched,
1554 default occurs when model_ids are passed instead,
1555 ie outside a CC modeling session
1556
1557 (default: [])
1558 @type star_grid: list[Star()]
1559 @keyword models: list of pacs_ids or gastronoom model ids, default if
1560 Star models are passed instead
1561
1562 (default: [])
1563 @type models: list[strings]
1564 @keyword exclude_data: if enabled only the sphinx mdels are plotted.
1565
1566 (default: 0)
1567 @type exclude_data: bool
1568 @keyword fn_plt: A plot filename to which an index is added for each
1569 subband.
1570
1571 (default: '')
1572 @type fn_plt: string
1573 @keyword cfg: path to the Plotting2.plotCols config file. If default, the
1574 hard-coded default plotting options are used.
1575
1576 (default: '')
1577 @type cfg: string
1578 @keyword fn_trans_marker: A file that includes TRANSITION definitions.
1579 These transitions will be marked up in the
1580 plot. For instance, when indicating a subset
1581 of transitions for one reason or another.
1582 The line type can be set for this specific
1583 subset, differently from other lines and
1584 regardless of the molecule. In combination
1585 with a doubly defined line label (through the
1586 star_grid['GAS_LINES']), lines can be marked
1587 up.
1588
1589 (default: '')
1590 @type fn_trans_marker: string
1591 @keyword mark_undetected: Mark the undetected transitions in the same
1592 way extra marked transitions would be marked
1593 by fn_trans_marker.
1594
1595 (default: 0)
1596 @type mark_undetected: bool
1597 @keyword include_band: Include a name tag for the band order in plot.
1598
1599 (default: 1)
1600 @type include_band: bool
1601 @keyword line_label_dashedvib: Use dashed lines for vibrational
1602 transitions in line label lines.
1603
1604 (default: 0)
1605 @type line_label_dashedvib: bool
1606
1607 '''
1608
1609 print '***********************************'
1610 print '** Creating PACS + Sphinx plot.'
1611 if self.pacs is None:
1612 print '** No PATH_PACS given. Cannot plot PACS spectra without data'+\
1613 ' information. Aborting...'
1614 return
1615 if not star_grid and models:
1616 star_grid = self.makeStars(models=models)
1617 elif (not models and not star_grid) or (models and star_grid):
1618 print '** Input is undefined or doubly defined. Aborting.'
1619 return
1620 if set([s['MOLECULE'] and 1 or 0 for s in star_grid]) == set([0]):
1621 return
1622
1623 self.setSphinxPacs(star_grid)
1624 print '** Plotting now...'
1625
1626 cfg_dict = Plotting2.readCfg(cfg)
1627 if cfg_dict.has_key('filename'):
1628 fn_plt = cfg_dict.pop('filename')
1629 if cfg_dict.has_key('exclude_data'):
1630 exclude_data = bool(cfg_dict['exclude_data'])
1631 if cfg_dict.has_key('fn_trans_marker'):
1632 fn_trans_marker = cfg_dict['fn_trans_marker']
1633 if cfg_dict.has_key('include_band'):
1634 include_band = bool(cfg_dict['include_band'])
1635 if cfg_dict.has_key('mark_undetected'):
1636 mark_undetected = cfg_dict['mark_undetected']
1637 if cfg_dict.has_key('labels'):
1638 labels = bool(cfg_dict['labels'])
1639 else:
1640 labels = []
1641
1642
1643 if mark_undetected:
1644 trl = Transition.extractTransFromStars(star_grid,dtype='PACS')
1645 for ifn in range(len(self.pacs.data_filenames)):
1646 self.pacs.intIntMatch(trans_list=trl,ifn=ifn)
1647
1648 lls = self.createLineLabels(star_grid=star_grid,\
1649 fn_trans_marker=fn_trans_marker,\
1650 mark_undetected=mark_undetected,\
1651 instrument='PACS')
1652
1653 plot_filenames = []
1654 for wave,flux,sphinx_flux,dfn,band in zip(self.pacs.data_wave_list,\
1655 self.pacs.data_flux_list,\
1656 self.sphinx_flux_list,\
1657 self.pacs.data_filenames,\
1658 self.pacs.data_ordernames):
1659
1660 pfn = fn_plt if fn_plt else os.path.split(dfn)[1]
1661 pfn = self.setFnPlt(pfn,fn_suffix=band,fn_subfolder='PACS_results')
1662
1663 keytags = ['Model %i: %s'%(i+1,str(star['LAST_PACS_MODEL'])\
1664 .replace('_','\_'))
1665 for i,star in enumerate(star_grid)]
1666 if exclude_data:
1667 x_list = [wave]*(len(sphinx_flux))
1668 y_list = sphinx_flux
1669 else:
1670 x_list = [wave]*(len(sphinx_flux)+1)
1671 y_list = [flux]+sphinx_flux
1672 keytags = ['PACS Spectrum'] + keytags
1673 if include_band:
1674 elabel = [(band,0.05,0.80)]
1675 else:
1676 elabel = []
1677
1678 plot_title = '{} - {}'.format(self.star_name_plots,band)
1679 plot_filenames.append(Plotting2.plotCols(x=x_list,y=y_list,\
1680 keytags=keytags,number_subplots=2,\
1681 plot_title=plot_title,\
1682 cfg=cfg_dict,\
1683 line_labels=lls,\
1684 histoplot=not exclude_data and [0] or [],\
1685 filename=pfn,labels=labels+elabel,\
1686 line_label_spectrum=1,line_label_color=1))
1687 if plot_filenames and plot_filenames[0][-4:] == '.pdf':
1688
1689 pfn = fn_plt if fn_plt else 'PACS_spectrum'
1690 newf = self.setFnPlt(pfn) + '.pdf'
1691 DataIO.joinPdf(old=sorted(plot_filenames),new=newf,\
1692 del_old=not fn_plt)
1693 print '** Your plots can be found at:'
1694 print newf
1695 print '***********************************'
1696 else:
1697 print '** Your plots can be found at:'
1698 print '\n'.join(plot_filenames)
1699 print '***********************************'
1700
1701
1702
1703 - def plotPacsSegments(self,star_grid,pacs_segments_path='',mode='sphinx',\
1704 fn_plt='',fn_trans_marker='',cfg='',\
1705 include_sphinx=None,exclude_data=0):
1706
1707 '''
1708 Plot segments of spectra only.
1709
1710 An inputfile gives the wavelength ranges, given by pacs_segments_path.
1711
1712 Can include the sphinx results overplotted with the data, as well as line
1713 labels generated either for sphinx results (mode == 'sphinx') or from a
1714 spectroscopic database (mode == 'll').
1715
1716 @param star_grid: star models for which PACS data will be fetched,
1717 @type star_grid: list(Star())
1718
1719 @keyword pacs_segments_path: The path to the file listing pairs of
1720 wavelength ranges for plotting the
1721 segments. This par can be passed through
1722 the cfg file as well.
1723 @type pacs_segments_path: string
1724 @keyword mode: the mode in which this method is used, the string is
1725 added to the outputfilename, can be 'sphinx' or 'll' for
1726 now, determines the type of line labels. 'll' gives line
1727 labels generated from a spectroscopic database. 'sphinx'
1728 gives line labels for all transitions in all Star
1729 objects in star_grid. Can be passed through the cfg file.
1730
1731 (default: sphinx)
1732 @type mode: string
1733 @keyword include_sphinx: Add sphinx results to the plots
1734
1735 (default: None)
1736 @type include_sphinx: bool
1737 @keyword exclude_data: The data are not included when True.
1738
1739 (default: 0)
1740 @type exclude_data: bool
1741 @keyword fn_plt: A plot filename to which an index is added for each
1742 subband.
1743
1744 (default: '')
1745 @type fn_plt: string
1746 @keyword fn_trans_marker: A file that includes TRANSITION definitions.
1747 These transitions will be marked up in the
1748 plot. For instance, when indicating a subset
1749 of transitions for one reason or another.
1750 The line type can be set for this specific
1751 subset, differently from other lines and
1752 regardless of the molecule. In combination
1753 with a doubly defined line label (through the
1754 star_grid['GAS_LINES']), lines can be marked
1755 up.
1756
1757 (default: '')
1758 @type fn_trans_marker: string
1759 @keyword cfg: path to the Plotting2.plotCols config file. If default,
1760 the hard-coded default plotting options are used.
1761
1762 (default: '')
1763 @type cfg: string
1764
1765 '''
1766
1767 cfg_dict = Plotting2.readCfg(cfg)
1768 if cfg_dict.has_key('mode'):
1769 mode = cfg_dict['mode']
1770 if cfg_dict.has_key('pacs_segments_path'):
1771 pacs_segments_path = cfg_dict['pacs_segments_path']
1772 if cfg_dict.has_key('include_sphinx'):
1773 include_sphinx = cfg_dict['include_sphinx']
1774 if cfg_dict.has_key('filename'):
1775 fn_plt = cfg_dict.pop('filename')
1776 if cfg_dict.has_key('fn_trans_marker'):
1777 fn_trans_marker = cfg_dict['fn_trans_marker']
1778 if cfg_dict.has_key('exclude_data'):
1779 exclude_data = bool(cfg_dict['exclude_data'])
1780
1781 if self.pacs is None:
1782 print 'No PACS data found for plotting PACS segments. Aborting...'
1783 return
1784 elif not pacs_segments_path:
1785 print 'No pacs_segments_path given. Pass in the cfg file or in ' +\
1786 'the method call. Aborting...'
1787 return
1788 else:
1789 self.setSphinxPacs(star_grid)
1790
1791 if mode == 'll':
1792 xmins=[min(wave_list) for wave_list in self.pacs.data_wave_list]
1793 xmaxs=[max(wave_list) for wave_list in self.pacs.data_wave_list]
1794 lls = self.createLineLabelsFromLineLists(star=star_grid[0],\
1795 xmin=min(xmins),\
1796 xmax=max(xmaxs),\
1797 fn_trans_marker=\
1798 fn_trans_marker,\
1799 instrument='PACS')
1800 elif mode == 'sphinx':
1801 lls = self.createLineLabels(star_grid=star_grid,\
1802 fn_trans_marker=fn_trans_marker,
1803 instrument='PACS')
1804 else:
1805 print 'Mode for plotting PACS segments not recognized. Aborting...'
1806 return
1807
1808 if include_sphinx is None:
1809 include_sphinx = self.sphinx_flux_list and 1 or 0
1810 print '** Plotting spectral segments.'
1811 for index, (wmin, wmax) in \
1812 enumerate(zip(*DataIO.readCols(pacs_segments_path))):
1813 delta = (wmax-wmin)/2.
1814 for i_file,(wave,flux,filename) in \
1815 enumerate(zip(self.pacs.data_wave_list,\
1816 self.pacs.data_flux_list,\
1817 self.pacs.data_filenames)):
1818 if wmin > wave[0] and wmax < wave[-1]:
1819 flux = flux[abs(wave-((wmax+wmin)/2.))<=delta]
1820 if not include_sphinx: sphinx_flux = []
1821 else:
1822 sphinx_flux = \
1823 [f[abs(wave-((wmax+wmin)/2.))<=delta]
1824 for f in self.sphinx_flux_list[i_file] if list(f)]
1825 wave = wave[abs(wave-((wmax+wmin)/2.))<=delta]
1826
1827
1828 pfn = fn_plt if fn_plt else 'PACS_spectrum'
1829 suff = '{}_segment_{:.2f}_{:.2f}'.format(mode,wmin,wmax)
1830 subf = mode=='ll' and 'LineLists' or 'PACS_results'
1831 pfn = self.setFnPlt(pfn,fn_suffix=suff,fn_subfolder=subf)
1832
1833
1834 extra_stats = dict([('line_labels',lls),\
1835 ('histoplot',not exclude_data \
1836 and [0] or []),\
1837 ('filename',pfn)])
1838
1839
1840 w = [wave]*(len(sphinx_flux)+(not exclude_data and 1 or 0))
1841 f = exclude_data and sphinx_flux or [flux]+sphinx_flux
1842
1843 pfn = Plotting2.plotCols(x=w,y=f,cfg=cfg_dict,**extra_stats)
1844 print '** Segment finished and saved at:'
1845 print pfn
1846
1847
1848
1849 - def plotSpire(self,star_grid=[],models=[],exclude_data=0,fn_plt='',\
1850 fn_trans_marker='',number_subplots=3,cfg=''):
1851
1852 '''
1853 Plot SPIRE data along with Sphinx results. In flux (Jy) vs wavelength
1854 (micron)
1855
1856 @keyword star_grid: star models for which SPIRE data will be fetched,
1857 default occurs when model_ids are passed instead,
1858 ie outside a CC modeling session
1859
1860 (default: [])
1861 @type star_grid: list[Star()]
1862 @keyword models: list of gastronoom model ids, default if
1863 Star models are passed instead
1864
1865 (default: [])
1866 @type models: list[strings]
1867 @keyword exclude_data: if enabled only the sphinx mdels are plotted.
1868
1869 (default: 0)
1870 @type exclude_data: bool
1871 @keyword fn_trans_marker: A file that includes TRANSITION definitions.
1872 These transitions will be marked up in the
1873 plot. For instance, when indicating a subset
1874 of transitions for one reason or another.
1875 The line type can be set for this specific
1876 subset, differently from other lines and
1877 regardless of the molecule. In combination
1878 with a doubly defined line label (through the
1879 star_grid['GAS_LINES']), lines can be marked
1880 up.
1881
1882 (default: '')
1883 @type fn_trans_marker: string
1884 @keyword fn_plt: A plot filename to which an index is added for each
1885 subband.
1886
1887 (default: '')
1888 @type fn_plt: string
1889 @keyword cfg: path to the Plotting2.plotCols config file. If default,
1890 the hard-coded default plotting options are used.
1891
1892 (default: '')
1893 @type cfg: string
1894
1895 '''
1896
1897 print '***********************************'
1898 print '** Creating SPIRE + Sphinx plot.'
1899 if self.spire is None:
1900 print '** No PATH_SPIRE given. Cannot plot SPIRE spectra '+\
1901 'without data information. Aborting...'
1902 return
1903 if not star_grid and models:
1904 star_grid = self.makeStars(models=models)
1905 elif (not models and not star_grid) or (models and star_grid):
1906 print '** Input is undefined or doubly defined. Aborting.'
1907 return
1908 if set([s['MOLECULE'] and 1 or 0 for s in star_grid]) == set([0]):
1909 return
1910 print '** Plotting now...'
1911
1912 cfg_dict = Plotting2.readCfg(cfg)
1913 if cfg_dict.has_key('filename'):
1914 fn_plt = cfg_dict.pop('filename')
1915 if cfg_dict.has_key('exclude_data'):
1916 exclude_data = bool(cfg_dict['exclude_data'])
1917 if cfg_dict.has_key('fn_trans_marker'):
1918 fn_trans_marker = cfg_dict['fn_trans_marker']
1919 self.spire.prepareSphinx(star_grid)
1920 lls = self.createLineLabels(star_grid,instrument='SPIRE')
1921
1922 if fn_trans_marker:
1923 used_indices = list(set([ll[-2] for ll in lls]))
1924 this_index = [ii for ii in range(100) if ii not in used_indices][0]
1925 ells = self.createLineLabels(fn_trans_marker=fn_trans_marker,\
1926 ilabel=this_index,instrument='SPIRE')
1927 lls = lls + ells
1928 plot_filenames = []
1929 for wave,flux,band,dfn in zip(self.spire.data_wave_list,\
1930 self.spire.data_flux_list,\
1931 self.spire.data_ordernames,\
1932 self.spire.data_filenames):
1933
1934 pfn = fn_plt if fn_plt else 'spire_spectrum'
1935 pfn = self.setFnPlt(pfn,fn_suffix=band)
1936
1937 sphinx_flux = [self.spire.getSphinxConvolution(star,dfn)[1]
1938 for star in star_grid]
1939 w = exclude_data \
1940 and [wave]*(len(sphinx_flux)) \
1941 or [wave]*(len(sphinx_flux)+1)
1942 f = exclude_data and sphinx_flux or [flux] + sphinx_flux
1943 keytags = ['Model %i: %s'%(i+1,star['LAST_GASTRONOOM_MODEL']\
1944 .replace('_','\_'))
1945 for i,star in enumerate(star_grid)]
1946 if not exclude_data:
1947 keytags = ['Spire Spectrum'] + keytags
1948 plot_filenames.append(Plotting2.plotCols(x=w,y=f,\
1949 keytags=keytags,number_subplots=3,\
1950 line_label_color=1,line_labels=lls,\
1951 plot_title='%s: %s' %(self.plot_id.replace('_','\_'),\
1952 self.star_name_plots),\
1953 histoplot= not exclude_data and [0] or [],\
1954 filename=pfn,cfg=cfg_dict,\
1955 line_label_spectrum=1))
1956 if plot_filenames and plot_filenames[0][-4:] == '.pdf':
1957
1958 pfn = fn_plt if fn_plt else 'SPIRE_spectrum'
1959 newf = self.setFnPlt(pfn) + '.pdf'
1960 DataIO.joinPdf(old=sorted(plot_filenames,reverse=True),new=newf,\
1961 del_old=not fn_plt)
1962 print '** Your plots can be found at:'
1963 print newf
1964 print '***********************************'
1965 else:
1966 print '** Your plots can be found at:'
1967 print '\n'.join(plot_filenames)
1968 print '***********************************'
1969
1970
1971 - def plotIntTmb(self,star_grid=[],scale=1,fn_plt='',cfg = ''):
1972
1973 '''
1974 Plot of the integrated main beam temperature of the molecular data
1975 and that obtained from model(s). Jup vs K.
1976
1977 @keyword star_grid: star models to be included in
1978
1979 (default: [])
1980 @type star_grid: list[Star()]
1981 @keyword scale: scale int. Tmb to an antenna of 1 m**2,
1982 necesarry to compare data from different telescope_string
1983
1984 (default = 1)
1985 @type scale: bool
1986 @keyword fn_plt: A plot filename to which an index is added for each
1987 subband.
1988
1989 (default: '')
1990 @type fn_plt: string
1991 @keyword cfg: path to the Plotting2.plotCols config file. If default,
1992 the hard-coded default plotting options are used.
1993
1994 (default: '')
1995 @type cfg: string
1996
1997 '''
1998
1999 print '***********************************'
2000 print '** Plotting integrated main beam temperatures'
2001 print '** Plots can be found at '
2002
2003
2004 cfg_dict = Plotting2.readCfg(cfg)
2005 if cfg_dict.has_key('filename'):
2006 fn_plt = cfg_dict.pop('filename')
2007 if cfg_dict.has_key('scale'):
2008 scale = bool(cfg_dict['scale'])
2009
2010
2011
2012
2013 trans = star_grid[0]['GAS_LINES']
2014 tele = list(set([t.telescope for t in trans]))
2015 molecs = list(set([t.molecule for t in trans]))
2016
2017 for jj in range(len(molecs)):
2018 jup_split = []
2019 data = []
2020 error = []
2021
2022 for ii in range(len(tele)):
2023 tr = [t for t in trans
2024 if t.telescope == tele[ii] and t.molecule == molecs[jj]]
2025 jup_split.append([t.jup for t in tr])
2026 idata,ierror = Transition.getLineStrengths(tr,mode='dtmb',\
2027 scale=scale)
2028 data.append(idata)
2029
2030 error.append(idata*ierror)
2031
2032 tra = [t for t in trans if t.molecule == molecs[jj]]
2033
2034
2035 types_data = ['ro','gs','bp']
2036 if len(jup_split)%3 != 0:
2037 types_data = types_data + types_data[0:(len(jup_split)%3)]
2038 if len(tele) == 2:
2039 types_data = ['ro','gs']
2040
2041
2042
2043 data_model = []
2044 label_model = []
2045 types_model = []
2046 colors = ['k', 'm', '0.50', 'c', 'y', 'r', 'g','b']
2047 C = len(colors)
2048 for ii in range(len(star_grid)):
2049 trans = star_grid[ii]['GAS_LINES']
2050 molecs = list(set([t.molecule for t in trans]))
2051 intra = [t for t in trans if t.molecule == molecs[jj]]
2052 mod = Transition.getLineStrengths(intra,mode='mtmb',scale=scale)
2053 data_model.append(mod[0])
2054 label_model.append(trans[0].getModelId().replace('_','\_'))
2055 types_model.append("".join([colors[ii%C], '--x']))
2056
2057
2058 jup = [t.jup for t in tra]
2059 indices = np.argsort(jup)
2060 jup.sort()
2061 data_model = [list(data_model[x][indices])
2062 for x in range(len(data_model))]
2063
2064
2065
2066
2067 x_toplot = jup_split + [jup]*len(star_grid)
2068 y_toplot = data + data_model
2069 yerr_toplot = error + [None]*len(star_grid)
2070 types = types_data + types_model
2071 labels = tele + label_model
2072
2073
2074 pfn = fn_plt if fn_plt else 'intTmb'
2075 suff = '_'.join([star_grid[0]['LAST_GASTRONOOM_MODEL'], \
2076 molecs[jj].makeLabel()[2:-2]])
2077 pfn = self.setFnPlt(pfn,fn_suffix=suff)
2078
2079
2080 extra_pars = dict()
2081 extra_pars['xmin'] = min(jup)-0.5
2082 extra_pars['xmax'] = max(jup)+0.5
2083 extra_pars['figsize'] = (15,9)
2084 extra_pars['yaxis'] = '$\int T_\mathrm{mb}\ (\mathrm{K\ km/s})$'
2085 extra_pars['xaxis'] = '$J_{up}$'
2086 extra_pars['keytags'] = labels
2087 extra_pars['line_types'] = types
2088 extra_pars['filename'] = pfn
2089
2090
2091 pfn = Plotting2.plotCols(x=x_toplot,y=y_toplot,cfg=cfg_dict,\
2092 yerr=yerr_toplot,**extra_pars)
2093 print pfn
2094 print '***********************************'
2095