# -*-coding:Utf-8 -*
# ====================================================================

# ====================================================================
#   Packages
# ====================================================================
import sys
import numpy as np
from contourmag.contourmag import contourmag
# ====================================================================
#   Functions
# ====================================================================
def magnifcalc(t, param, Ds=None, tb=None):
    """Return the quadrupolar approximation of the magnification."""
### Get parameters
    t0 = params['t0']
    u0 = params['u0']
    tE = params['tE']
    rho = params['rho']
    gamma = params['gamma']
    q = params['q']
    piEN = params['piEN']
    piEE = params['piEE']
    alpha0 = params['alpha']
    s0 = params['s']
    dalpha = params['dalpha']
    ds = params['ds']
### Lens orbital motion
    alpha, s = lens_rotation(alpha0, s0, dalpha, ds, t, tb)
### Parallax
    DsN = Ds['N']
    DsE = Ds['E']
    tau = (t-t0)/tE + piEN * DsN + piEE * DsE
    beta = u0 + piEN * DsE - piEE * DsN
    x, y = binrot(alpha, tau, beta)
### Conversion secondary body left -> right
    x = -x
### Compute magnification
    accuracy = 1.e-2 # Absolute mag accuracy (mag+/-accuracy)
    return contourmag(s, q, rho, gamma, x, y, accuracy)
# --------------------------------------------------------------------
def binrot(theta, x_old, y_old):
    """Rotation by an angle alpha.

    :param theta: float, angle in radians.
    :param x_old: numpy array, x coodinate in the old frame.
    :param y_old: numpy array, y coodinate in the old frame.
    :return x_new: numpy array, x coodinate in the new frame.
    :return y_new: numpy array, y coodinate in the new frame.
    """
    cos_theta = np.cos(theta)
    sin_theta = np.sin(theta)
    x_new = x_old * cos_theta - y_old * sin_theta
    y_new = x_old * sin_theta + y_old * cos_theta
    return x_new, y_new
# --------------------------------------------------------------------
def lens_rotation(alpha0, s0, dalpha, ds, t, tb):
    """Compute the angle alpha and projected separation s for each
    time step due to the lens orbital motion.

    :param alpha0: angle alpha at date tb.
    :param s0: projected separation at date tb.
    :param dalpha: float, angular velocity at date tb
        (radians.year^-1).
    :param ds: change rate of separation (year^-1).
    :param t: list of dates.
    :param tb: time reference for linear development.
    :type alpha0: float
    :type s0: float
    :type dalpha: float
    :type ds: float
    :type t: numpy array
    :type tb: float
    :return: unpacked list of actual alpha and s values at each date.
    :rtype: numpy array, numpy array
    """
    Cte_yr_d = 365.25  # Julian year in days
    alpha = alpha0 - (t - tb) * dalpha / Cte_yr_d
    s = s0 + (t-tb) * ds / Cte_yr_d
    return alpha, s


@@ -271,16 +271,75 @@
 y as np%0A
+from muLAn.models.contourmag.contourmag import contourmag%0A#
 from con
@@ -375,16 +375,16 @@
 tourmag%0A
-
 # ======
@@ -677,17 +677,16 @@
  = param
-s
 %5B't0'%5D%0A 
@@ -698,17 +698,16 @@
  = param
-s
 %5B'u0'%5D%0A 
@@ -719,17 +719,16 @@
  = param
-s
 %5B'tE'%5D%0A 
@@ -741,17 +741,16 @@
  = param
-s
 %5B'rho'%5D%0A
@@ -766,17 +766,16 @@
  = param
-s
 %5B'gamma'
@@ -789,17 +789,16 @@
  = param
-s
 %5B'q'%5D%0A  
@@ -807,25 +807,24 @@
 piEN = param
-s
 %5B'piEN'%5D%0A   
@@ -836,17 +836,16 @@
  = param
-s
 %5B'piEE'%5D
@@ -863,17 +863,16 @@
  = param
-s
 %5B'alpha'
@@ -891,14 +891,41 @@
 aram
-s
 %5B's'%5D
+ # no orbital motion s0 -%3E s
 %0A   
@@ -939,17 +939,16 @@
  = param
-s
 %5B'dalpha
@@ -964,17 +964,16 @@
  = param
-s
 %5B'ds'%5D%0A#
@@ -1380,16 +1380,17 @@
  return 
+%5B
 contourm
@@ -1393,16 +1393,19 @@
 ourmag(s
+%5Bi%5D
 , q, rho
@@ -1414,19 +1414,25 @@
 gamma, x
-, y
+%5Bi%5D, y%5Bi%5D
 , accura
@@ -1434,16 +1434,40 @@
 ccuracy)
+ for i in range(len(x))%5D
 %0A# -----

Update ESBLcont.py