Python pylab.average函数代码示例

def flow_rate_hist(sheets):
    ant_rates = []
    weights = []
    for sheet in sheets:
        ants, seconds, weight = flow_rate(sheet)
        ant_rate = seconds / ants
        #ant_rate = ants / seconds
        ant_rates.append(ant_rate)
        weights.append(float(weight))
        #weights.append(seconds)

    weights = pylab.array(weights)
    weights /= sum(weights)

    #print "ants per second"
    print "seconds per ant"
    mu = pylab.mean(ant_rates)
    print "mean", pylab.mean(ant_rates)
    wmean = pylab.average(ant_rates, weights=weights)
    print "weighted mean", wmean
    print "median", pylab.median(ant_rates)
    print "std", pylab.std(ant_rates, ddof=1)
    ant_rates = pylab.array(ant_rates)
    werror = (ant_rates - mu) * weights
    print "weighted std", ((sum(werror ** 2))) ** 0.5
    print "weighted std 2", (pylab.average((ant_rates - mu)**2, weights=weights)) ** 0.5
    pylab.figure()
    pylab.hist(ant_rates)
    pylab.savefig('ant_flow_rates.pdf', format='pdf')
    pylab.close()

def makeplot(filename):
    T0 = 2452525.374416
    P = 0.154525
    
    X = pl.load(filename)
    x = X[:,0]
    y = X[:,1]
    print x[0] # check for HJD faults
    
    #orbital phase
    p = (x-T0)/P
    
    pl.figure(figsize=(6,4))
    pl.subplots_adjust(hspace=0.47,left=0.16)
    
    pl.subplot(211)
    pl.scatter(p,y,marker='o',s=0.1,color='k')
    pl.ylim(-0.06,0.06)
    pl.xlim(pl.average(p)-1.25,pl.average(p)+1.25)
    pl.ylabel('Intensity')
    pl.xlabel('Orbital Phase')
    
    pl.subplot(212)
    f,a = ast.signal.dft(x,y,0,4000,1)
    pl.plot(f,a,'k')
    pl.ylabel('Amplitude')
    pl.xlabel('Frequency (c/d)')
    #pl.ylim(yl[0],yl[1])
    
    #pl.vlines(3636,0.002,0.0025,color='k',linestyle='solid')
    #pl.vlines(829,0.002,0.0025,color='k',linestyle='solid')
    #pl.text(3500,0.00255,'DNO',fontsize=11)
    #pl.text(700,0.00255,'lpDNO',fontsize=11)
    pl.ylim(0.0,0.004)
    pl.savefig('%spng'%filename[:-3])

def AT(s_ij=0, s_ik=0, s_jk=0, S=0):
	'''
	Calculates: SEFD = (2k/S) * (s_ij*s_ik)/(s_jk)
	'''
	kb = 1.38e3 # Boltzmann's Constant in Jy m^2 K^-1
	s_ij = pl.average(s_ij)
	s_ik = pl.average(s_ik)
	s_jk = pl.average(s_jk)
	return (2*kb/S)*(s_ij*s_ik)/(s_jk-s_ij*s_ik)

def visualize ():
    sample_rate, snd = load_sample(".\\hh-closed\\dh9.WAV")
    print snd.dtype
    data = normalize(snd)
    print data.shape
    n = data.shape[0]
    length = float(n)
    print length / sample_rate, "s"
    timeArray = arange(0, length, 1)
    timeArray = timeArray / sample_rate
    timeArray = timeArray * 1000  #scale to milliseconds
    ion()
    if False:
        plot(timeArray, data, color='k')
        ylabel('Amplitude')
        xlabel('Time (ms)')
        raw_input("press enter")
        exit()
    p = fft(data) # take the fourier transform
    nUniquePts = ceil((n+1)/2.0)
    print nUniquePts
    p = p[0:nUniquePts]
    p = abs(p)
    p = p / float(n) # scale by the number of points so that
                 # the magnitude does not depend on the length
                 # of the signal or on its sampling frequency
    p = p**2  # square it to get the power

    # multiply by two (see technical document for details)
    # odd nfft excludes Nyquist point
    if n % 2 > 0: # we've got odd number of points fft
        p[1:len(p)] = p[1:len(p)] * 2
    else:
        p[1:len(p) -1] = p[1:len(p) - 1] * 2 # we've got even number of points fft

    print p
    freqArray = arange(0, nUniquePts, 1.0) * (sample_rate / n);
    plot(freqArray/1000, 10*log10(p), color='k')
    xlabel('Frequency (kHz)')
    ylabel('Power (dB)')
    raw_input("press enter")

    m = average(freqArray, weights = p)
    v = average((freqArray - m)**2, weights= p)
    r = sqrt(mean(data**2))
    s = var(data**2)
    print "mean freq", m #TODO: IMPORTANT: this is currently the mean *power*, not the mean freq.  What we want is mean freq weighted by power
    print "var freq", v
    print "rms", r
    print "squared variance", s

def show_grey_channels(I):
    K = average(I, axis=2)
    for i in range(3):
        J = zeros_like(I)
        J[:, :, i] = K
        figure(i+10)
        imshow(J)

def hist_values(parameter, group, strategy, decay_type, label, y_limit=None):
    values = group[parameter]
    values = list(values)
    
    binsize = 0.05
    if 'zoom' in label:
        binsize = 0.01
    
    if y_limit == None:
        y_limit = max(values)
    cutoff = 1
    #weights = np.ones_like(values)/float(len(values))
    
    weights = group['num_lines'] / sum(group['num_lines'])
    weights = np.array(weights)
    
    mu = pylab.average(values, weights=weights)
    sigma2 = pylab.var(values)
    
    pylab.figure()
    pylab.hist(values, weights=weights, bins=np.arange(0, cutoff + binsize, binsize))
    title_items = []
    title_items.append('%s maximum likelihood values %s %s %s' % (parameter, strategy, decay_type, label))
    title_items.append('mean of estimates = %f' % mu)
    title_items.append('variance of estimates = %f' % sigma2)
    title_str = '\n'.join(title_items)
    #pylab.title(parameter + ' maximum likelihood values ' + str(strategy) + ' ' + str(outname))
    #pylab.title(title_str)
    print title_str
    pylab.xlabel('%s mle' % parameter, fontsize=20)
    pylab.ylabel('weighted proportion', fontsize=20)
    pylab.xlim((0, 1))
    pylab.ylim((0, y_limit))
    pylab.savefig('repair_ml_hist_%s_%s_%s_%s.pdf' % (parameter, strategy, decay_type, label), format='pdf')
    pylab.close()

 def movingaverage(x,L):
     ma = pl.zeros(len(x),dtype='Float64')
     # must take the lead-up zone into account (prob slow)
     for i in range(0,L):
         ma[i] = pl.average(x[0:i+1])
 
     for i in range(L,len(x)):
         ma[i] = ma[i-1] + 1.0/L*(x[i]-x[i-L])
         
     return ma

def img2ascii(filename, map_array=None):
    a = imread(filename)

    print "Converting ..."
    # useful only when reading .jpg files.
    # PIL is used for jpegs; converting PIL image to numpy array messes up. 
    # a = a[::-1, :] 

    # convert image to grayscale.
    if len(a.shape) > 2:
        a = 0.21 * a[:,:,0] + 0.71 * a[:,:,1] + 0.07 * a[:,:,2]
    a_r, a_c = a.shape[:2]
    a_max = float(a.max())

    blk_siz = 1 #size of block

    if map_array == None:
        # just linearly map gray level to characters.
        # giving lowest gray level to space character.
        out_file = open(filename + 'lin' + str(blk_siz) + '.txt', 'w')
        print "File %s opened" %out_file.name
        for i in range(0, a_r, blk_siz*2):
            for j in range(0, a_c, blk_siz):
                b = a[i:i+2*blk_siz, j:j+blk_siz]
                b_char = chr(32+int((1-average(b))*94))
                out_file.write(b_char)
            out_file.write("\n")
        out_file.close()
    else:
        # map based on visual density of characters.
        out_file = open(filename + 'arr' + str(blk_siz) + '.txt', 'w')
        print "File %s opened" %out_file.name
        for i in range(0, a_r, blk_siz*2):
            for j in range(0, a_c, blk_siz):
                b = a[i:i+2*blk_siz, j:j+blk_siz]
                b_mean = int(average(b)/a_max*(len(map_array)-1))
                b_char = chr(map_array[b_mean])
                out_file.write(b_char)
            out_file.write("\n")
        out_file.close()
        
    print "%s Converted! \nWritten to %s" %(filename, out_file.name)

def sigclip(im,nsig):
    # returns min and max values of image inside nsig sigmas
    temp = im.ravel()
    sd = pl.std(temp)
    m = pl.average(temp)
    gt = temp > m-nsig*sd
    lt = temp < m+nsig*sd
    temp = temp[gt*lt]
    mini = min(temp)
    maxi = max(temp)
    
    return mini,maxi

    def plot(self):
        """generate the plot formatting"""
        if self.data == None:
            print "Must load and parse data first"
            sys.exit()
            
        for k,v in self.data.iteritems():
            for type, data in v.iteritems():
                pylab.clf()
                height = int(self.height)
                width = int(self.width)
                pylab.figure()
                ax = pylab.gca()
                ax.set_xlabel('<--- Width = %s wells --->' % str(width))
                ax.set_ylabel('<--- Height = %s wells --->' % str(height))
                ax.set_yticks([0,height/10])
                ax.set_xticks([0,width/10])
                ax.set_yticklabels([0,height])
                ax.set_xticklabels([0,width])
                ax.autoscale_view()
                pylab.jet()
            #color = self.makeColorMap()
            #remove zeros for calculation of average
                flattened = []
                for i in data:
                    for j in i:
                        flattened.append(j)
                flattened = filter(lambda x: x != 0.0, flattened)
                Aver=pylab.average(flattened)
                name = type.replace(" ", "_")
                fave = ("%.2f") % Aver
                pylab.title(k.strip().split(" ")[-1] + " Heat Map (Average = "+fave+"%)")
                ticks = None
                vmax = None
                if type == "Region DR":
                    ticks = [0.0,0.2,0.4,0.6,0.8,1.0]
                    vmax = 1.0
                else:
                    ticks = [0.0,0.4,0.8,1.2,1.6,2.0]
                    vmax = 2.0
                    
                pylab.imshow(data, vmin=0, vmax=vmax, origin='lower')
                pylab.colorbar(format='%.2f %%',ticks=ticks)
                pylab.vmin = 0.0
                pylab.vmax = 2.0
            #pylab.colorbar()     

                if self.savePath is None:
                    save = "%s_heat_map_%s.png" % (name,k)
                else:
                    save = path.join(self.savePath,"%s_heat_map_%s.png" % (name,k))
                pylab.savefig(save)
                pylab.clf()

def process_window(sample_rate, data):
    # print "processing window"
    # print data.dtype
    # print data.shape
    n = data.shape[0]
    length = float(n)
    # print length / sample_rate, "s"
    p = fft(data) # take the fourier transform
    nUniquePts = ceil((n+1)/2.0)
    p = p[0:nUniquePts]
    p = abs(p)
    p = p / float(n) # scale by the number of points so that
                 # the magnitude does not depend on the length
                 # of the signal or on its sampling frequency
    p = p**2  # square it to get the power

    # multiply by two (see technical document for details)
    # odd nfft excludes Nyquist point
    if n % 2 > 0: # we've got odd number of points fft
        p[1:len(p)] = p[1:len(p)] * 2
    else:
        p[1:len(p) -1] = p[1:len(p) - 1] * 2 # we've got even number of points fft
    freqArray = arange(0, nUniquePts, 1.0) * (sample_rate / n);

    if sum(p) == 0:
        raise Silence
    m = average(freqArray, weights = p)
    v = sqrt(average((freqArray - m)**2, weights= p))
    r = sqrt(mean(data**2))
    s = var(data**2)
    print "mean freq", m #TODO: IMPORTANT: this is currently the mean *power*, not the mean freq.  What we want is mean freq weighted by power
    # print freqArray
    # print (freqArray - m)
    # print p
    print "var freq", v
    print "rms", r
    print "squared variance", s
    return [m, v, r, s]

 def writeMetricsFile(self, filename):
     '''Writes the lib_cafie.txt file'''
     if self.data == None:
         print "Must load and parse data first"
         sys.exit()
     cafie_out = open(filename,'w')
     for k,v in self.data.iteritems():
         for type, data in v.iteritems():
             flattened = []
             for i in data:
                 for j in i:
                     flattened.append(j)
             flattened = filter(lambda x: x != 0.0, flattened)
             Aver=pylab.average(flattened)
             name = type.replace(" ", "_")
             if 'LIB' in k:
                 if len(flattened)==0:
                     cafie_out.write('%s = %s\n' % (name,0.0))
                     Aver = 0
                 else:
                     cafie_out.write('%s = %s\n' % (name,Aver))
                     Aver=pylab.average(flattened)
     cafie_out.close()

 def calcTDData(self,tdDatas):
     #tdDatas is a a 3d array of measurements, along with their uncertainties
     #meantdData is the weighted sum of the different measurements
     #meantdData,sumofweights=py.average(tdDatas[:,:,1:3],axis=0,weights=1.0/tdDatas[:,:,3:]**2,returned=True)
     meantdData=py.average(tdDatas[:,:,1:3],axis=0)
     #use error propagation formula
     noise=py.sqrt(py.mean(self.getAllPrecNoise()[0]**2))
     if tdDatas.shape[0]==1:
         rep = py.zeros((len(tdDatas[0,:,0]),2))
     else:
         rep = py.std(tdDatas[:,:,1:3],axis=0, ddof=1)/py.sqrt(self.numberOfDataSets)
     unc = py.sqrt(rep**2+noise**2)
     #unc=py.sqrt(1.0/sumofweights)
     #time axis are all equal
     return py.column_stack((tdDatas[0][:,0],meantdData,unc))       

def amptime(uv, baseline="0_1", pol="xx", applycal=False):
	'''
	Plots Amp vs Time for a single baseline. 
	'''
	fig = pl.figure()
	ax = fig.add_subplot(111)
	aipy.scripting.uv_selector(uv, baseline, pol)
	for preamble, data, flags in uv.all(raw=True):
		uvw, t, (i, j) = preamble
		ax.plot(t, pl.average(pl.absolute(data)), 'ks', mec='None', alpha=0.2, ms=5)
	hfmt = dates.DateFormatter('%m/%d %H:%M')
	ax.xaxis.set_major_locator(dates.HourLocator())
	ax.xaxis.set_major_formatter(hfmt)
	ax.set_ylim(bottom = 0)
	pl.xticks(rotation='vertical')	

def calculateScores(arr, HEIGHT, WIDTH):
    rowlen,collen = arr.shape
    scores = pylab.zeros(((rowlen/INCREMENT),(collen/INCREMENT)))
    score = []
    for row in range(rowlen/INCREMENT):
        for column in range(collen/INCREMENT):            
            keypassed,size = getAreaScore(row,column,arr, HEIGHT, WIDTH)
            scores[row,column] = round(float(keypassed)/float(size)*100,2)
            if keypassed > 2:
                score.append(round(float(keypassed)/float(size)*100,2))         
            
            #scores[0,0] = 0
            #scores[HEIGHT/INCREMENT -1,WIDTH/INCREMENT -1] = 100
    print scores
    
    flattened = []
    for i in score:
        flattened.append(i)
    flattened = filter(lambda x: x != 0.0, flattened)    
    average=pylab.average(flattened)
    
    return score, scores, average