def int_pars(x, y, vds, **kwargs):
    """
     calculates York regression and Coe parameters (with Tauxe Fvds)
    """
    # first do linear regression a la York
    # do Data Model 3 way:
    if 'version' in list(kwargs.keys()) and kwargs['version'] == 3:
        n_key = 'int_n_measurements'
        b_key = 'int_b'
        sigma_key = 'int_b_sigma'
        f_key = 'int_f'
        fvds_key = 'int_fvds'
        g_key = 'int_g'
        q_key = 'int_q'
        b_beta_key = 'int_b_beta'

    else:  # version 2
        n_key = 'specimen_int_n'
        b_key = 'specimen_b'
        sigma_key = 'specimen_b_sigma'
        f_key = 'specimen_f'
        fvds_key = 'specimen_fvds'
        g_key = 'specimen_g'
        q_key = 'specimen_q'
        b_beta_key = 'specimen_b_beta'

    xx, yer, xer, xyer, yy, xsum, ysum, xy = 0., 0., 0., 0., 0., 0., 0., 0.
    xprime, yprime = [], []
    pars = {}
    pars[n_key] = len(x)
    n = float(len(x))
    if n <= 2:
        print("shouldn't be here at all!")
        return pars, 1
    for i in range(len(x)):
        xx += x[i]**2.
        yy += y[i]**2.
        xy += x[i] * y[i]
        xsum += x[i]
        ysum += y[i]
    xsig = np.sqrt(old_div((xx - (old_div(xsum**2., n))), (n - 1.)))
    ysig = np.sqrt(old_div((yy - (old_div(ysum**2., n))), (n - 1.)))
    sum = 0
    for i in range(int(n)):
        yer += (y[i] - old_div(ysum, n))**2.
        xer += (x[i] - old_div(xsum, n))**2.
        xyer += (y[i] - old_div(ysum, n)) * (x[i] - old_div(xsum, n))
    slop = -np.sqrt(old_div(yer, xer))
    pars[b_key] = slop
    s1 = 2. * yer - 2. * slop * xyer
    s2 = (n - 2.) * xer
    sigma = np.sqrt(old_div(s1, s2))
    pars[sigma_key] = sigma
    s = old_div((xy - (xsum * ysum / n)), (xx - old_div((xsum**2.), n)))
    r = old_div((s * xsig), ysig)
    pars["specimen_rsc"] = r**2.
    ytot = abs(old_div(ysum, n) - slop * xsum / n)
    for i in range(int(n)):
        xprime.append(old_div((slop * x[i] + y[i] - ytot), (2. * slop)))
        yprime.append((old_div((slop * x[i] + y[i] - ytot), 2.)) + ytot)
    sumdy, dy = 0, []
    dyt = abs(yprime[0] - yprime[int(n) - 1])
    for i in range((int(n) - 1)):
        dy.append(abs(yprime[i + 1] - yprime[i]))
        sumdy += dy[i]**2.
    f = old_div(dyt, ytot)
    pars[f_key] = f
    pars["specimen_ytot"] = ytot
    ff = old_div(dyt, vds)
    pars[fvds_key] = ff
    ddy = (old_div(1., dyt)) * sumdy
    g = 1. - old_div(ddy, dyt)
    pars[g_key] = g
    q = abs(slop) * f * g / sigma
    pars[q_key] = q
    pars[b_beta_key] = old_div(-sigma, slop)
    return pars, 0