Skip to content

Instantly share code, notes, and snippets.

@chrisgorgo

chrisgorgo/AFNI_ANTs_affine_transforms.ipynb

Created September 14, 2018 16:06
Show Gist options
  • Save chrisgorgo/ffcab0cb0fa3298e2196e0d5c7d6cd3c to your computer and use it in GitHub Desktop.
Save chrisgorgo/ffcab0cb0fa3298e2196e0d5c7d6cd3c to your computer and use it in GitHub Desktop.
Download ZIP
Display the source blob
Display the rendered blob
Raw
AFNI_ANTs_affine_transforms.ipynb
{
"cells": [
{
"cell_type": "code",
"execution_count": 136,
"metadata": {},
"outputs": [],
"source": [
"from scipy.io import loadmat, savemat\n",
"from numpy import loadtxt, around, hstack, vstack, zeros, float64"
]
},
{
"cell_type": "code",
"execution_count": 86,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{'AffineTransform_double_3_3': array([[ 9.99771194e-01],\n",
" [-4.22685205e-05],\n",
" [-1.66161657e-05],\n",
" [-4.15907894e-06],\n",
" [ 8.65957305e-01],\n",
" [-4.99806282e-01],\n",
" [-7.38437417e-05],\n",
" [ 4.99559786e-01],\n",
" [ 8.65365078e-01],\n",
" [-1.17959414e-03],\n",
" [ 8.09175497e+00],\n",
" [ 4.96808506e+00]]), 'fixed': array([[0.],\n",
" [0.],\n",
" [0.]])}"
]
},
"execution_count": 86,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"loadmat('sub-001_desc-deobliquedrotatedrealignedants_T1w0GenericAffine.mat')"
]
},
{
"cell_type": "code",
"execution_count": 56,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{'AffineTransform_double_3_3': array([[ 9.99771194e-01],\n",
" [-4.22685205e-05],\n",
" [-1.66161657e-05],\n",
" [-4.15907894e-06],\n",
" [ 8.65957305e-01],\n",
" [-4.99806282e-01],\n",
" [-7.38437417e-05],\n",
" [ 4.99559786e-01],\n",
" [ 8.65365078e-01],\n",
" [-1.17959414e-03],\n",
" [ 8.09175497e+00],\n",
" [ 4.96808506e+00]]), 'fixed': array([[0.],\n",
" [0.],\n",
" [0.]])}"
]
},
"execution_count": 56,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"loadmat('sub-001_desc-deobliquedrotatedrealignedants_T1w0GenericAffine.mat')"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [],
"source": [
"orig_ants_mat = loadmat('sub-001_desc-deobliquedrotatedrealignedants_T1w0GenericAffine.mat')['AffineTransform_double_3_3']"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [],
"source": [
"orig_afni_mat = loadtxt('sub-001_desc-deobliquedrotatedrealigned_T1w.aff12.1D')"
]
},
{
"cell_type": "code",
"execution_count": 22,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 1. , -0. , -0. , -0. ],\n",
" [ 0.87, -0.5 , -0. , 0.5 ],\n",
" [ 0.87, -0. , 8.09, 4.97]])"
]
},
"execution_count": 22,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"around(orig_ants_mat, 2).reshape(3,4, order='C')"
]
},
{
"cell_type": "code",
"execution_count": 23,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 1. , -0. , -0. , -0. ],\n",
" [-0. , 0.87, 0.5 , 8.09],\n",
" [-0. , -0.5 , 0.87, 4.97]])"
]
},
"execution_count": 23,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"around(orig_ants_mat, 2).reshape(3,4, order='F')"
]
},
{
"cell_type": "code",
"execution_count": 24,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 1. , -0. , -0. , -0. ],\n",
" [ 0.87, -0.5 , -0. , 0.5 ],\n",
" [ 0.87, -0. , 8.09, 4.97]])"
]
},
"execution_count": 24,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"around(orig_ants_mat, 2).reshape(3,4, order='A')"
]
},
{
"cell_type": "code",
"execution_count": 21,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 1. , -0. , 0. , 0. ],\n",
" [-0. , 0.87, -0.5 , 8.08],\n",
" [ 0. , 0.5 , 0.87, 4.96]])"
]
},
"execution_count": 21,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"around(orig_afni_mat, 2).reshape(3,4, order='C')"
]
},
{
"cell_type": "code",
"execution_count": 25,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 1. , 0. , -0.5 , 0.5 ],\n",
" [-0. , -0. , 8.08, 0.87],\n",
" [ 0. , 0.87, 0. , 4.96]])"
]
},
"execution_count": 25,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"around(orig_afni_mat, 2).reshape(3,4, order='F')"
]
},
{
"cell_type": "code",
"execution_count": 31,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 1. , -0. , 0. , 0. ],\n",
" [-0. , 0.87, -0.5 , 8.08],\n",
" [ 0. , 0.5 , 0.87, 4.96]])"
]
},
"execution_count": 31,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"around(orig_afni_mat, 2).reshape(3,4, order='C')"
]
},
{
"cell_type": "code",
"execution_count": 35,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 1. , -0. , -0. , -0. ],\n",
" [ 0.87, -0.5 , -0. , 0.5 ],\n",
" [ 0.87, -0. , 8.09, 4.97]])"
]
},
"execution_count": 35,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"around(orig_ants_mat, 2).reshape(3,4, order='C')"
]
},
{
"cell_type": "code",
"execution_count": 34,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 1. , -0. , -0. ],\n",
" [-0. , 0.87, -0.5 ],\n",
" [-0. , 0.5 , 0.87]])"
]
},
"execution_count": 34,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"around(orig_ants_mat[:9], 2).reshape(3,3, order='C')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"```\n",
"mat=[reshape(afftransform(1:9),[3,3])',afftransform(10:12)];\n",
"\n",
"m_Translation=mat(:,4);\n",
"mat=[mat;[0,0,0,1]];\n",
"\n",
"\n",
"for i=1:3\n",
" m_Offset(i) = m_Translation(i) + m_Center(i);\n",
" for j=1:3\n",
" m_Offset(i) = m_Offset(i)-(mat(i,j) * m_Center(j)); % (i,j) should remain the same since in C indexing rows before cols, too.\n",
" end\n",
"end\n",
"\n",
"mat(1:3,4)=m_Offset;\n",
"mat=inv(mat);\n",
"\n",
"% convert RAS to LPS (ITK uses LPS)\n",
"mat=mat.*...\n",
" [1 1 -1 -1\n",
" 1 1 -1 -1\n",
" -1 -1 1 1\n",
" 1 1 1 1];\n",
"```"
]
},
{
"cell_type": "code",
"execution_count": 140,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 1. , -0. , -0. , -0. ],\n",
" [-0. , 0.87, -0.5 , 8.09],\n",
" [-0. , 0.5 , 0.87, 4.97],\n",
" [ 0. , 0. , 0. , 1. ]])"
]
},
"execution_count": 140,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def read_ants_affine(input_file, debug=False):\n",
" \n",
" def calculateOffset(affine_transform, translations, fixed):\n",
" offset = zeros(3)\n",
" for i in range(3):\n",
" offset[i] = translation[i] + fixed[i]\n",
" for j in range(3):\n",
" offset[i] = offset[i] - (affine_transform[i,j] * fixed[j])\n",
" return offset.reshape(3,-1)\n",
"\n",
" mat_content = loadmat(input_file)\n",
" if debug:\n",
" print(mat_content)\n",
" matching_affine_keys = [k for k in mat_content.keys() if k.startswith('AffineTransform')]\n",
" assert len(matching_affine_keys) == 1\n",
" matching_affine_key = matching_affine_keys[0]\n",
" affine_transform = mat_content[matching_affine_key][:9].reshape(3,3, order='C')\n",
" translation = mat_content[matching_affine_key][9:]\n",
" fixed = mat_content['fixed']\n",
" offset = calculateOffset(affine_transform, translations, fixed)\n",
" combined = vstack((hstack((affine_transform, offset)),[0,0,0,1]))\n",
" return combined\n",
" \n",
"around(read_ants_affine('sub-001_desc-deobliquedrotatedrealignedants_T1w0GenericAffine.mat'), 2)"
]
},
{
"cell_type": "code",
"execution_count": 144,
"metadata": {},
"outputs": [],
"source": [
"def read_afni_affine(input_file, debug = False):\n",
" orig_afni_mat = loadtxt(input_file)\n",
" if debug:\n",
" print(orig_afni_mat)\n",
"\n",
" combined = vstack((orig_afni_mat.reshape(3,4, order='C'), [0,0,0,1]))\n",
" return combined"
]
},
{
"cell_type": "code",
"execution_count": 89,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([0. , 8.08, 4.96])"
]
},
"execution_count": 89,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"around(read_afni_affine('sub-001_desc-deobliquedrotatedrealigned_T1w.aff12.1D'), 2)[:3,3]"
]
},
{
"cell_type": "code",
"execution_count": 74,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 9.99771194e-01, -4.22685205e-05, -1.66161657e-05,\n",
" -1.17959414e-03],\n",
" [-4.15907894e-06, 8.65957305e-01, -4.99806282e-01,\n",
" 8.09175497e+00],\n",
" [-7.38437417e-05, 4.99559786e-01, 8.65365078e-01,\n",
" 4.96808506e+00]])"
]
},
"execution_count": 74,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"hstack((affine_transform, calculateOffset(affine_transform, translations, fixed)))"
]
},
{
"cell_type": "code",
"execution_count": 141,
"metadata": {},
"outputs": [],
"source": [
"def write_ants_affine(affine, out_file, debug=False):\n",
" out_dict = {}\n",
" out_dict['AffineTransform_double_3_3'] = hstack((affine[:3,:3].reshape(1, -1), affine[:3,3].reshape(1, -1))).reshape(-1,1).astype(float64)\n",
" out_dict['fixed'] = zeros((3,1))\n",
" if debug:\n",
" print(out_dict)\n",
" \n",
" savemat(out_file, out_dict, format='4')"
]
},
{
"cell_type": "code",
"execution_count": 153,
"metadata": {
"scrolled": true
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"[ 0.999991 -0.00110987 -0.00419122 0.0665791 0.00107965 0.999973\n",
" -0.00720623 -0.0969192 0.00419911 0.00720163 0.999965 1.06055 ]\n",
"{'AffineTransform_double_3_3': array([[ 0.999991 ],\n",
" [-0.00110987],\n",
" [-0.00419122],\n",
" [ 0.00107965],\n",
" [ 0.999973 ],\n",
" [-0.00720623],\n",
" [ 0.00419911],\n",
" [ 0.00720163],\n",
" [ 0.999965 ],\n",
" [ 0.0665791 ],\n",
" [-0.0969192 ],\n",
" [ 1.06055 ]]), 'fixed': array([[0.],\n",
" [0.],\n",
" [0.]])}\n"
]
}
],
"source": [
"write_ants_affine(read_afni_affine('C:/stanford_gdrive/data/transformation_test_data/vol0001tovol0020.aff12.1D', True), \n",
" \"test.mat\", True)"
]
},
{
"cell_type": "code",
"execution_count": 132,
"metadata": {},
"outputs": [],
"source": [
"savemat?"
]
},
{
"cell_type": "code",
"execution_count": 155,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"(238, 12)"
]
},
"execution_count": 155,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"loadtxt(\"C:/stanford_gdrive/data/transformation_test_data/volred.aff12.1D\").shape"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.6.6"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment