# Least square fitting between spac coefficient with bessel function zero order first kind

I used to applying least squares fitting of spac coefficient with Bessel function zero order first kind such that according to book by Hiroshi Okada, The Microtremor Survey Method (2003). I used scipy.optimize.least_square. I want to implement equation

Here is my script:

    # DEFINING DATA
r = 135
x =[]
y=[]
for i,j in zip(avgspaccoef['Frequency (Hz)'], avgspaccoef['SPAC Coefficient']):
if 0.8< i <= 6:
x.append(i)
y.append(j)

# RESAMPLE DATA TO PREVENT VERY DENSE ABCYSSA (IT REFERS TO FREQUENCY)
x_interp = np.linspace(np.round(min(x), decimals=1), np.round(max(x), decimals=1), 100)
y_interp=np.interp(x_interp, x,y)
x=x_interp
y=y_interp

# PERFORM POLYNOMIAL APROXIMATION TO SMOOTH THE COARSE SPAC COEFFICIENT DATA
coefficient_polyaprx= np.polyfit(x, y, 13)
result= np.polyval(coefficient_polyaprx, x)
y_fit = []
for i in x:
y_value = 0
for coeff, pow in zip(coefficient_polyaprx, range(len(coefficient_polyaprx)-1, -1, -1)):
y_value += coeff * i**pow
y_fit.append(y_value)

# IMPLEMENTING LEAST SQUARE TO SOLVE XI
def res(xi, fi, r, y_fit):
omega=2*np.pi*r
return y_fit-jv(0,(omega*fi/xi))

def calculate_dispcurv(x, y, r, res, x0):
pv = []
f_pv = []
bessel_fit = []
x_bessel = []
for i in x:
fi = i
index = np.argmin(np.abs(np.array(x) - fi))
fi = x[index]
yi = y[index]

#Using Least Square
fit_y = least_squares(res, x0, args=(fi, r, yi), method='dogbox', bounds=(100,1200), loss='linear')
result = np.abs(fit_y.x)

pv.append(result)
f_pv.append(fi)
bessel_fit.append(yi - fit_y.fun)
x_bessel.append(2 * np.pi * r * fi / result)
return pv, f_pv, bessel_fit, x_bessel

x0=[200, 300, 500, 600, 700, 2000]
pv, f_pv, bessel_fit, x_bessel = calculate_dispcurv(x_interp, y_interp, r ,res, x0)


The script will yield results similiar to following:

*note that data and result used in this example was gained from Asten MW, Hayashi K. 2018. Application of the Spatial Auto-Correlation Method for Shear-Wave Velocity Studies Using Ambient Noise. Surveys in Geophysics [Internet]. 39(4):633–659. https://doi.org/10.1007/s10712-018-9474-2

My script tends to perfoming each index in x separately and perform the least squares computation. In order to yield robust results it is proper to compute least squares between frequency and spac coefficient as data and Bessel function where xi is unknown variable. Does anyone have ideas on how to help me figure it out, please?

• Is this more relevant to post on a python or computation SE forum? Commented 11 hours ago