import os
import sys
# Get the directory containing the current file
current_dir = os.path.dirname(os.path.abspath(__file__))
sys.path.append(os.path.join(current_dir, '../buildBG/'))
from buildBG import *
from utilities import *
import math
import numpy as np

file_path=os.path.join(current_dir, '../Electrogenic cotransporter/')
fmatrix=file_path+'SLC5_f.csv'
rmatrix=file_path+'SLC5_r.csv'
bg_components_json='../buildBG/BG_components.json'        

# build the bond graph model
bg_dict={}       
direction = 'e2f'
add_BG_components(bg_dict,bg_components_json, fmatrix)
add_BG_connections(bg_dict, fmatrix, direction)
direction = 'f2e'
add_BG_connections(bg_dict, rmatrix, direction)

# update the equations of the bond graph model
voi={'description': 'Time', 'units': 'second', 'symbol': 't'}
update_BG_eqn(bg_dict,voi)

fmatrix=file_path+'SLC5_f_chem.csv'
rmatrix=file_path+'SLC5_r_chem.csv' 
# update the BG parameters for the biochemical reactions
eName, eID, ePort, fName, fID, fPort,N_f=load_matrix(fmatrix)
eName, eID, ePort, fName, fID, fPort,N_r=load_matrix(rmatrix)
AVO=6.022e23
V1=1e-3

# Parameters for Fig 5, k_61=3

k_12=8e4*V1*V1
k_23=1e5*V1
k_34=50
k_45=800
k_56=10
k_61=3
k_25=0.3
k_21=500
k_32=20
k_43=50
k_65=50*V1*V1
k_16=35
k_52=k_12*k_25*k_56*k_61/(k_21*k_65*k_16)
k_54=k_23*k_34*k_45*k_52/(k_32*k_43*k_25)
balance1=k_12*k_23*k_34*k_45*k_56*k_61/(k_21*k_32*k_43*k_54*k_65*k_16)
balance2=k_12*k_25*k_56*k_61/(k_21*k_52*k_65*k_16)
balance3=k_23*k_34*k_45*k_52/(k_32*k_43*k_54*k_25)
print(balance1,balance2,balance3)
kf=[k_12, k_23, k_34, k_45, k_56, k_61, k_25]
kr=[k_21, k_32, k_43, k_54, k_65, k_16, k_52 ]
K_c=[]
N_c=[]
V_E=1
V_o=0.085
V_i=0.085
Ws=[V_i,V_o,V_i,V_o,V_E,V_E,V_E,V_E,V_E,V_E]
kappa, K, K_eq, diff_,  zero_est,k_est= kinetic2BGparams(N_f,N_r,kf,kr,K_c,N_c,Ws)
# print the estimated results
for i in range(len(kappa)):
    print(fName[i],kappa[i])
for i in range(len(K)):
    print(eName[i],K[i])
k_est_name=['k_12','k_23','k_34','k_45','k_56','k_61','k_25','k_21','k_32','k_43','k_54','k_65','k_16','k_52']
for i in range(len(k_est)):
    print(k_est_name[i],k_est[i])
print('sum of relative errors',diff_)
update_BioBG_params(bg_dict, kappa, fName, K, eName)
update_BG_params(bg_dict,[('T',293)])
# combine K and kappa to a single array
param_val=np.concatenate((kappa,K))
param_fName=[]
for i in fName:
    param_fName.append('$\kappa_{'+i+'}$')
param_eName=[]
for i in eName:
    param_eName.append('$K_{'+i+'}$')
param_name=param_fName+param_eName
param_units=['fmol$.s^{-1}$']*len(fName)+['fmol$^{-1}$']*len(eName)
write_params_csv(param_name,param_val,param_units,csv_file=file_path+'SLC5_BG.csv')
# save the bond graph model to a json file   
json_file=file_path+'SLC5_BG.json'    
save_json(bg_dict, json_file) 

# Parameters for Fig 10, k_61=5

# Force the thermodynamic constants to be the same as the Fig 5
K_Na=3.216e-1
K_Glc=4.848e1
K_1=2.235
K_2=10.437
K_3=8.602

k_12=8e4*V1*V1
k_23=1e5*V1
k_34=50
k_45=800
k_56=10
k_61=5
k_25=0.3
k_21=500
k_32=20
k_43=50
k_65=50*V1*V1
k_16=35
k_52=k_12*k_25*k_56*k_61/(k_21*k_65*k_16)
k_54=k_23*k_34*k_45*k_52/(k_32*k_43*k_25)
balance1=k_12*k_23*k_34*k_45*k_56*k_61/(k_21*k_32*k_43*k_54*k_65*k_16)
balance2=k_12*k_25*k_56*k_61/(k_21*k_52*k_65*k_16)
balance3=k_23*k_34*k_45*k_52/(k_32*k_43*k_54*k_25)
print(balance1,balance2,balance3)
kf=[k_12, k_23, k_34, k_45, k_56, k_61, k_25]
kr=[k_21, k_32, k_43, k_54, k_65, k_16, k_52 ]
V_E=1
V_o=0.085
V_i=0.085
Ws=[V_i,V_o,V_i,V_o,V_E,V_E,V_E,V_E,V_E,V_E]
K_c=[K_Na*V_i,K_Glc*V_i,K_1,K_2,K_3]
N_c=[[1,0,0,0,0,0,0,0,0,0],[0,0,1,0,0,0,0,0,0,0],[0,0,0,0,1,0,0,0,0,0],[0,0,0,0,0,1,0,0,0,0],[0,0,0,0,0,0,1,0,0,0]]
kappa, K, K_eq, diff_,  zero_est,k_est= kinetic2BGparams(N_f,N_r,kf,kr,K_c,N_c,Ws)
# print the estimated results
for i in range(len(kappa)):
    print(fName[i],kappa[i])
for i in range(len(K)):
    print(eName[i],K[i])
k_est_name=['k_12','k_23','k_34','k_45','k_56','k_61','k_25','k_21','k_32','k_43','k_54','k_65','k_16','k_52']
for i in range(len(k_est_name)):
    print(k_est_name[i],k_est[i])
print('sum of relative errors',diff_)
update_BioBG_params(bg_dict, kappa, fName, K, eName)
update_BG_params(bg_dict,[('T',293)])
# combine K and kappa to a single array
param_val=np.concatenate((kappa,K))
param_fName=[]
for i in fName:
    param_fName.append('$\kappa_{'+i+'}$')
param_eName=[]
for i in eName:
    param_eName.append('$K_{'+i+'}$')
param_name=param_fName+param_eName
param_units=['fmol$.s^{-1}$']*len(fName)+['fmol$^{-1}$']*len(eName)
write_params_csv(param_name,param_val,param_units,csv_file=file_path+'SLC5_BG_Fig10.csv')
# save the bond graph model to a json file   
json_file=file_path+'SLC5_BG_Fig10.json'    
save_json(bg_dict, json_file)