Monday, April 4, 2016

A simple multi-panel velocity plot

def velocity_plot(wave, fnu, line_label, line_center, redshift, vwin, Ncol) :
    '''Make a plot of flux density versus rest-frame velocity for several spectral lines.
    Inputs are:
    wave:        observed wavelength (currently Angstrom)
    fnu:         flux density (erg/s/cm^2/Hz)  ** should this be flambda instead?
    line_label:  tuple of line labels
    line_center: np array of rest-frame wavelengths (Angstroms)
    redshift:    redshift to assume,
    vwin:        velocity window to use (km/s)
    Ncol:        number of columns to plot.  Number of rows set automatically.
    '''
    Nrow = ceil(line_center.size / Ncol)  # Calculate how many rows to generate
    fig = plt.figure()
    restwave = wave / (1.0 + redshift)

    for ii, dum in enumerate(line_label) :
        print "Trying to plot ", line_label[ii], " at ", line_center[ii]
        ax = fig.add_subplot(Nrow, Ncol, ii+1)
        vel =  (restwave - line_center[ii])/line_center[ii] * A_c/1E5     # km/s   
        in_window = np.where( (vel >= vwin*-1 ) & (vel <= vwin) )
        plt.step(vel[in_window], fnu[in_window])       
        plt.plot( (-1*vwin, vwin), (1.0,1.0))  # plot unity continuum
        plt.ylim(0.0, 1.3)  # May need to change these limits
        plt.xlim(-1*vwin, vwin)
        if ii == len(line_label) -1 :
            plt.xlabel("rest-frame velocity (km/s)")  # if last subplot, make xlabel
        else :
            ax.axes.xaxis.set_ticklabels([])  # if not last subplot, suppress  numbers on x axis
        plt.annotate( line_label[ii], (0.3,0.8), xycoords="axes fraction")
#        plt.ylabel("relative fnu")
        fig.subplots_adjust(hspace=0)

No comments:

Post a Comment