Python numpy.ln函数代码示例

def pdf_x(L, l, delta, dv):
    if -L / 2 + l / 2 + dv < delta < L / 2 - l / 2 - dv:
        return 2 * dv / l * ln(L / (L - l))
    if  L / 2 - l / 2 - dv < delta < L / 2 - l / 2 + dv:
        return 1 / 2. + 1 / l * ((delta - dv) * ln(L - l) + (delta + dv) * (1 - ln(2 * (delta + dv))) + 2 * dv * ln(L) - L / 2.)
    if (delta > L / 2. - l / 2. + dv) and (delta > L / 2. - l / 2. - dv):
        return 1 / l * (2 * dv * (1 + ln(L / 2.)) - (delta + dv) * ln(delta + dv) + (delta - dv) * ln(delta - dv))

def _get_I_total(R0, d, t):
    A = np.ln(R0)
    B = np.ln(1+ d) 
    mu = 0.5 * A / B
    I = 0.5 * np.exp( 0.25 * A**2 / B * np.sqrt( np.pi / B))
    I = I * np.sqrt(B) * (np.erf(t - mu)- np.erf(-mu))
    return I

def GKB_1(u_zref, zref, h, LAI, Wfol, Ta, pa):
   """Same as FKB_1, but then for spatial in- and output"""
   
   # Constants
   C_d = 0.2   # foliage drag coefficient
   C_t = 0.05  # heat transfer coefficient
   k = 0.41     # Von Karman constant
   Pr = 0.7    # Prandtl number
   hs = 0.009  # height of soil roughness obstacles (0.009-0.024)
   
   # Calculations
   Wsoil = 1.0 - Wfol
   h = ifthenelse(Wfol == 0.0, hs, h)
   
   z0 = 0.136 * h   # Brutsaert (1982)
   u_h0 = u_zref * ln(2.446) / ln ((zref - 0.667 * h)/z0) # wind speed at canopy height
   ust2u_h = 0.32 - 0.264 / exp(15.1 * C_d * LAI)
   ustarh = ust2u_h * u_h0
   nu0 = 1.327E-5 * (101325.0/pa) * (Ta / 273.15 + 1.0) ** 1.81 # kinematic viscosity
   n_h = C_d * LAI / (2.0 * ust2u_h ** 2.0)
   # First term
   F1st = ifthenelse(pcrne(n_h, 0.0), k * C_d / (4.0 * C_t * ust2u_h * (1.0 - exp(pcrumin(n_h)/2.0))) * Wfol ** 2.0, 0.0)
   # Second term
   S2nd = k * ust2u_h * 0.136 * Pr ** (2.0/3.0) * sqrt(ustarh * h / nu0) * Wfol ** 2.0 * Wsoil ** 2.0
   # Third term
   T3rd = (2.46 * (u_zref * k / ln(zref/hs) * hs / nu0) ** 0.25 - ln(7.4)) * Wsoil ** 2.0
   
   return F1st + S2nd + T3rd

def pdf_x_0(L, l, delta, dv):
    if -L / 2 + l / 2 + dv < delta < L / 2 - l / 2 - dv:
        return 1 / l * ln(L / (L - l))
    if  L / 2 - l / 2 - dv < delta < L / 2 - l / 2 + dv:
        return 0
    if (delta > L / 2. - l / 2. + dv) and (delta > L / 2. - l / 2. - dv):
        return 1 / l * ln(L / (2 * delta))

def probability_cut_nooverlaps(lc, lf, delta):
    """
        Probability that we cut fiber of nooverlapped fibers
    """
    # type I and III
    if (lf <= lc / 2.0 - delta) and (lf <= lc / 2.0 + delta):
        # print "Varianta I or III ",
        return lc / lf * (ln(lc / (lc - lf))) - 1
    # type II
    if (lf < lc / 2.0 + delta) and (lf > lc / 2.0 - delta):
        # print "Varianta II ",
        return (
            1
            / 2.0
            / lf
            * (
                (lc + 2.0 * delta) * ln(2.0 / (lc + 2.0 * delta))
                - (ln(lc - lf) + 1) * (lc - 2.0 * delta)
                + 2.0 * lc * ln(lc)
            )
        )
    # type IV
    if (lf >= lc / 2.0 + delta) and (lf >= lc / 2.0 - delta):
        # print "Varianta IV ",
        return 1 + 1 / lf * (
            -lc + (lc / 2.0 - delta) * ln((lc + 2 * delta) / (lc - 2 * delta)) + lc * ln(2 * lc / (lc + 2 * delta))
        )

def loadData(path="../data/",k=5,log='add',pca_n=0,SEED=34):
	from pandas import DataFrame, read_csv
	from numpy import log as ln
	from sklearn.cross_validation import KFold
	from sklearn.preprocessing import LabelEncoder
	from sklearn.preprocessing import StandardScaler
	train = read_csv(path+"train.csv")
	test = read_csv(path+"test.csv")
	id = test.id
	target = train.target
	encoder = LabelEncoder()
	target_nnet = encoder.fit_transform(target).astype('int32')
	feat_names = [x for x in train.columns if x.startswith('feat')]
	train = train[feat_names].astype(float)
	test = test[feat_names]
	if log == 'add':
		for v in train.columns:
			train[v+'_log'] = ln(train[v]+1)
			test[v+'_log'] = ln(test[v]+1)
	elif log == 'replace':
		for v in train.columns:
			train[v] = ln(train[v]+1)
			test[v] = ln(test[v]+1)      
	if pca_n > 0:
		from sklearn.decomposition import PCA
		pca = PCA(pca_n)
		train = pca.fit_transform(train)
		test = pca.transform(test)
	scaler = StandardScaler()
	scaler.fit(train)
	train = DataFrame(scaler.transform(train),columns=['feat_'+str(x) for x in range(train.shape[1])])
	test = DataFrame(scaler.transform(test),columns=['feat_'+str(x) for x in range(train.shape[1])])
	cv = KFold(len(train), n_folds=k, shuffle=True, random_state=SEED)
	return train, test, target, target_nnet, id, cv, encoder

def Ar(z, lf):
    #z is the dist between the 
    lf = np.float(lf)
    #######################
            #int#
    #######################
    
    if type(z) == int:
        if z == 0:
            return 1
        if z == lf / 2:
            return 0
        return lf * (lf - 2. * z + 2. * z * (.6931471806 + ln((1. / lf) * z))) / (lf - 2. * z) ** 2
    
    #######################
            #ARRAY#
    #######################
    
    res = 2 * z / lf * (ln(2 * z / lf) - 1) + 1
    if any(z == lf / 2):
        z = list(z)
        ind_h = z.index(lf / 2)
        res[ind_h] = 0
        z = np.array(z)
    if any(z == 0):
        z = list(z)
        ind_z = z.index(0) 
        res[ind_z] = 1 
        z = np.array(z)
    return res 

def PSIma(f, g):
   a = 0.33
   b = 0.41
   pi = 3.141592654
   tangens = scalar(atan((2.0 * g - 1.0) / sqrt(3.0))) * pi /180
   tangens = ifthenelse(tangens > pi/2.0, tangens - 2.0 * pi, tangens) 
   PSIma = ln(a + f) - 3.0 * b * f ** (1.0 / 3.0) + b * a ** (1.0 / 3.0) / 2.0 * ln((1 + g) ** 2.0 / (1.0 - g + sqrt(g))) + sqrt(3.0) * b * a ** (1.0 / 3.0) * tangens
   return PSIma

def Cw(hi, L, z0, z0h):
   alfa = 0.12
   beta = 125.0
   C0 = (alfa / beta) * hi
   C1 = pcrumin(z0h) / L
   C11 = -alfa * hi / L
   C21 = hi / (beta * z0)
   C22 = -beta * z0 / L
   C = ifthenelse(z0 < C0, pcrumin(ln(alfa)) + PSIh_y(C11) - PSIh_y(C1), ln(C21) + PSIh_y(C22) - PSIh_y(C1))
   Cw = ifthenelse(C < 0.0, 0.0, C) # This results from unfortunate parameter combination!
   return Cw

def le(L, l):
    '''
        Solve embedded length l_e of fibes (including null values) (integral)
    '''
    if L < l:
        #print 'very short specimen',
        return L / 4.
    if L < 2. * l:
        #print 'short specimen',
        return -L / 4. - L / 4. * ln(L / 2.) + L / 4. * ln(L) + 1 / 2. * l * (1 - L / (2. * l)) + 1 / 4. * L * ln(L / (2. * l)) + l / 4.
    if L >= 2. * l:
        #print 'long specimen',
        return -1 / 4. * l - 1 / 4. * L * ln(L - l) + 1 / 4. * L * ln(L)

def prob( lc, lf, delta ):
    #type I and III
    if ( lf <= lc / 2. - delta ) and ( lf <= lc / 2. + delta ):
        print "Varianta I or III ",
        return lc / lf * ( ln( lc / ( lc - lf ) ) ) - 1
    #type II
    if ( lf < lc / 2. + delta ) and ( lf > lc / 2. - delta ):
        print "Varianta II ",
        return 1 / 2. / lf * ( ( lc + 2. * delta ) * ln( 2. / ( lc + 2. * delta ) ) - ( ln( lc - lf ) + 1 ) * ( lc - 2. * delta ) + 2. * lc * ln( lc ) )
    #type IV
    if ( lf >= lc / 2. + delta ) and ( lf >= lc / 2. - delta ):
        print "Varianta IV ",
        return  1 + 1 / lf * ( -lc + ( lc / 2. - delta ) * ln( ( lc + 2 * delta ) / ( lc - 2 * delta ) ) + lc * ln( 2 * lc / ( lc + 2 * delta ) ) )

def u_pbl(NDVI):
   """Calculates Planetary Boundary Layer wind speed [m s-1] from NDVI
   
   NDVI Input PCRaster NDVI map (scalar, ratio between 0 and 1)"""
   
   z0m = 0.005 + 0.5 * (nd_mid/nd_max) ** 2.5
   assert z0m >= 0.0
   fc = ((nd_mid - nd_min) / nd_df) ** 2.0    # fractional vegetation cover == Wfol (-)
   assert fc >= 0.0
   h = z0m / 0.136                            # total height of vegetation (m)
   d = 2.0/3.0 * h			      # zero plane displacement (m)
   u_c = ln((z_pbl - d) / z0m) / ln((z_ms - d) / z0m)
   u_pbl = u_s * u_c
   return u_pbl, z0m, d, fc, h

 def _get_results(self):
     """dividing the interval <a,b>,
     returns aprox x, error estimation, No. of steps ..."""
     int = [self.a, self.b]
     if self.f(self.a) * self.f(self.b) > 0:
         print "None or more than 1 roots in selected interval"
     else:
         while abs(int[0] - int[1])/2 > self.error:
             if self.f(int[0]) * self.f((int[0] + int[1])*0.5) < 0:
                 int.insert(1,(int[0]+int[1])*0.5), int.pop()
             else:
                 int.insert(1,(int[0]+int[1])*0.5), int.pop(0) 
         return [(int[0] + int[1])/2,
                 abs(int[0] - int[1])/2,
                 (ln(2)-ln( (abs( int[0] - int[1] ) /2)/(self.b-self.a)))/ln(2)]

def calc_bethe_entropy(B):
    """Calculate the Bethe entropy given beliefs B"""
    Sbethe = 0
    for roti,rotj in rotedges:
        try:
            if B[(roti,rotj)]>0:
                Sbethe -= B[(roti,rotj)]* ln(B[(roti,rotj)])
        except RuntimeError:
            pass  # can't find edge. . .
    sumqi = sum([len(res2partners[resid])-1 for resid in resids])
    blogb = 0
    for roti in eg.GetVertexIDs():
        if B[roti]>0:
            blogb += B[roti]* ln(B[roti])
    Sbethe += sumqi*blogb
    return Sbethe

def calc_meanfield_entropy(B):
    """Calculate the mean-field entropy given beliefs B"""
    Smeanfield = 0
    for roti in eg.GetVertexIDs():
        if B[roti]>0:
            Smeanfield -= B[roti] * ln(B[roti])
    return Smeanfield