#!/usr/bin/env python
'''
-------------------
Copyright (c) 2016 Computational Biomechanics (CoBi) Core, Department of
Biomedical Engineering, Cleveland Clinic
Permission is hereby granted, free of charge, to any person obtaining a
copy of this software and associated documentation files (the
"Software"), to deal in the Software without restriction, including
without limitation the rights to use, copy, modify, merge, publish,
distribute, sublicense, and/or sell copies of the Software, and to permit
persons to whom the Software is furnished to do so, subject to the
following conditions:
The above copyright notice and this permission notice shall be included
in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN
NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,
DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR
OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE
USE OR OTHER DEALINGS IN THE SOFTWARE.
-------------------
modify_model.py
Written by Craig Bennetts re-using scripts by Scott Sibole.
Updated and customized by Ahmet Erdemir.
Contact:
Ahmet Erdemir, PhD
Computational Biomodeling (CoBi) Core
Department of Biomedical Engineering
Lerner Research Institute
Cleveland Clinic
Cleveland, OH 44195, USA
erdemira@ccf.org
'''
import sys
import string
import xml.dom.minidom as xmd
import numpy as np
def usage(argv):
print
print 'USAGE: ' + argv[0] + ' <model.feb> <modify_model.cfg> <modified_model_filename.feb>'
print
def main(argv):
# Check for appropriate usage
if len(argv) != 4:
usage(argv)
sys.exit(1)
model_filename = argv[1]
cfg_filename = argv[2]
# set output filename
output_filename = argv[3]
# set clinical translations zero
clinical_dof = {'AP': 0.0,
'ML': 0.0,
'DC': 0.0,
'FE': 0.0,
'IE': 0.0,
'VV': 0.0}
# read config file
cfg_content = parse_file(cfg_filename, '*', '**')
# set clinical translations to values from config file
for key in cfg_content.keys():
if key in clinical_dof.keys():
clinical_dof[key] = float(cfg_content[key][0])
# calculate Grood & Suntay rotations and translations (motion of tibia relative to femur)
grood_dof = cli2grd(clinical_dof)
# calculate FEA rotations and translations (motion of femur relative to tibia)
fea_dof = grd2fea(grood_dof)
# read model file
model = xmd.parse(model_filename)
# write modified model file with new FEA rotations and translations
model_elements = model.documentElement
model_kinematics = model_elements.getElementsByTagName("prescribed")
for kmt in model_kinematics:
kmt.childNodes[0].data = fea_dof[kmt.getAttribute("bc")]
fout = open(output_filename,"wb")
model.writexml(fout, encoding=model.encoding)
fout.close()
# END OF main()
def parse_file(fname, skyw, scmt): # file name, keyword string, comment string
fid = open(fname, 'r')
# remove trailing white space and load file content
flines = map(string.rstrip, fid.readlines()) # each line is a string in a list
fid.close()
#remove spaces
tmp = []
dmys = flines[:]
for line in dmys:
dmy = line.replace(" ","")
tmp.append(dmy)
del flines
flines = tmp # same list with white spaces removed
del dmys
# remove comments and empty lines
dmys = flines[:]
for dmy in dmys:
if (scmt in dmy) or (dmy == ''):
flines.remove(dmy) # list with all comment lines and empty lines removed
del dmys
# find all keywords
findx = [];
i = 0;
for fline in flines:
if skyw in fline:
findx.append(i) # list of the line #s that contain keywords
i = i + 1
# save contents to a dictionary
fcontent = {}
nfind = len(findx) # get number of keywords
for i in range(nfind):
indst = findx[i]+1 # line number just after keyword
if i < nfind-1:
indsp = findx[i+1]; # line number of next keyword
else:
indsp = len(flines)
fcontent[flines[indst-1].strip(skyw)] = flines[indst:indsp] # remove * from keyword, use as key in dictionary, following list (up to next keyword) is the value
return fcontent
# END OF parse_file()
def cli2grd(clinical_dof):
# motion of tibia relative to femur
# based on clinical dof convention and signs
# based on model coordinate systems
grood_dof = {}
grood_dof['Rx'] = -clinical_dof['FE']
grood_dof['Ry'] = clinical_dof['VV']
grood_dof['Rz'] = clinical_dof['IE']
sRy = np.sin(np.deg2rad(grood_dof['Ry']))
grood_dof['x'] = (clinical_dof['ML'] - sRy*clinical_dof['DC'])/(-1.0 + sRy**2)
grood_dof['y'] = clinical_dof['AP']
grood_dof['z'] = (sRy*clinical_dof['ML'] - clinical_dof['DC'])/(1.0 - sRy**2)
return grood_dof
def grd2fea(grood_dof):
# motion of femur relative to tibia
# based on inverse of grood & suntay convention and signs
# to prescribe motion of femur for finite element analysis
# based on model coordinate systems
Tft = groodTransformation(grood_dof)
Ttf = np.linalg.inv(Tft)
theta = np.arccos((np.trace(Ttf[0:-1,0:-1])-1.)/2.)
if np.abs(theta) > 0:
s = 1./(2.*np.sin(theta))*np.array([Ttf[2][1] - Ttf[1][2],Ttf[0][2]-Ttf[2][0],Ttf[1][0]-Ttf[0][1]])
else:
s = np.array([0., 0., 0.])
fea_dof = {}
fea_dof['x'] = Ttf[0,3]
fea_dof['y'] = Ttf[1,3]
fea_dof['z'] = Ttf[2,3]
fea_dof['Rx'] = theta*s[0]
fea_dof['Ry'] = theta*s[1]
fea_dof['Rz'] = theta*s[2]
return fea_dof
def groodTransformation(grood_dof):
# motion of femur relative to tibia
# based on model coordinate systems
x = grood_dof['x']
y = grood_dof['y']
z = grood_dof['z']
Rx = np.deg2rad(grood_dof['Rx'])
Ry = np.deg2rad(grood_dof['Ry'])
Rz = np.deg2rad(grood_dof['Rz'])
Tx = np.array([[1., 0., 0., x], [0., np.cos(Rx), -np.sin(Rx), 0.], [0., np.sin(Rx), np.cos(Rx), 0.], [0.,0.,0.,1.]])
Ty = np.array([[np.cos(Ry), 0., np.sin(Ry), 0.],[0., 1., 0., y], [-np.sin(Ry), 0., np.cos(Ry), 0.], [0.,0.,0.,1.]])
Tz = np.array([[np.cos(Rz), -np.sin(Rz), 0., 0.], [np.sin(Rz), np.cos(Rz), 0., 0.],[0., 0., 1., z], [0.,0.,0.,1.]])
return np.dot(Tx, np.dot(Ty,Tz))
if __name__ == "__main__":
main(sys.argv)