Introduction About Site Map

XML
RSS 2 Feed RSS 2 Feed
Navigation

Main Page | Blog Index

Sunday, November 21st, 2010, 12:01 pm

Grids and Tagging – Octave Code Samples

Objects collection - stethoscope

FOR a little bit of background on tagging and my current work, see my previous posts on the subject [1, 2, 3], prior to me getting data from Professor Axel. The goal is to advance existing methods that include tagging by applying to them an algorithm which I developed for non-rigid registration assessment in the brain. This time the data that I deal with is cardiac, so temporal sets are a far more complicated problem, particularly due to persistent motion. Of interest there’s the paper “Single Heartbeat Cardiac Tagging for the Evaluation of Transient Phenomena” from Daniel A. Herzka, J. Andrew Derbyshire, Peter Kellman, and Elliot R. McVeigh. It can be found in Magnetic Resonance in Medicine (Volume 54, Issue 6, pages 1455–1464, December 2005) and the full paper says in its abstract: “Many cardiac abnormalities are of a transient nature, creating a beat-to-beat variation in myocardial function. This work presents the cardiac imaging technique for the measurement of regional function during transient cardiac phenomena. All information necessary for the reconstruction of a cine loop is acquired within a single heartbeat, avoiding the temporal blurring introduced by segmented imaging due to the assumption of cardiac cycle periodicity…”

I shall assume that motion can be treated like changes that occur over time on different dates or across subjects. The part in the above paper which deals with tags is: “Images reconstructed from such data displayed tag blurring and reduced tag persistence due to motion and interheartbeat variability. Images acquired during the Valsalva maneuver demonstrated apparent beat-to-beat variability, visible both in the images and as changing strain patterns and ventricular volumes.”

Therein lies a real problem because the tags cannot be assured to exist where they are needed the most. Fortunately, the literature contains many examples of work that can reconstruct imaged objects based on spatial information other than the tags themselves. A paper which was cited a lot in the 90s comes from Young and Axel (published in 1992, but date of current version is 6th of August 2002) and it is titled “Non-rigid heart wall motion using MR tagging”. Prior to it, in 1988, E A Zerhouni, D M Parish, W J Rogers, A Yang, and E P Shapiro looked at ways of analysing myocardial motion, potentially in troubled patients. Their paper “Human heart: tagging with MR imaging–a method for noninvasive assessment of myocardial motion” is summarised thusly:

Specified regions of the myocardium can be labeled in magnetic resonance (MR) imaging to serve as markers during contraction. The technique is based on locally perturbing the magnetization of the myocardium with selective radio-frequency (RF) saturation of multiple, thin tag planes during diastole followed by conventional, orthogonal-plane imaging during systole. The technique was implemented on a 0.38-T imager and tested on phantoms and volunteers. In humans, tags could be seen 60-450 msec after RF saturation, thus permitting sampling of the entire contractile phase of the cardiac cycle. Tagged regions appear as hypointense stripes, and their patterns of displacement reflect intervening cardiac motion. In addition to simple translation and rotation, complex motions such as cardiac twist can be demonstrated. The effects of RF pulse angle, relaxation times, and heart rate on depiction of the tagged region are discussed.

The German Heart Institute Berlin published the paper “Magnetic resonance stress tagging in ischemic heart disease” many years later (2005) and they used a 1.5-Tesla scanner for probing patients with coronary artery disease (CAD) under high-dose dobutamine stress. This clinical study helped show how “quantitative myocardial tagging may become a tool that reduces the need for high-dose dobutamine stress.”

As stated in previous posts and short papers which will be uploaded at a later date, my intention is to consider stress/strain analysis applied to the imaged tags. It is the tracking of the tags and making use of them that’s very challenging around boundaries that move throughout the cardiac cycle. Identifying templates like a cross where the tags have an intersection is a computationally expensive process that may also require taking uncertainty into account. If the goal is to track tags in real-time-level speeds (like in some ultrasound packages), then a faster approach is needed, e.g. shuffle distance. This approach would involve comparing different frames (at least one prior frame, with the exception of the first which has no predecessors) and then, in the junctions of the tags, best matches will be found and recorded in a vector of {x,y} coordinates, where these coordinates are refined over time and help quickly determine the motion of objects inside the image. Calculations using mechanical models (based upon shapes that the tags form) are an area which is already well explored, so it is the tracking and recovery of tags which my code deals with first and foremost. prototyping is done in GNU Octave, which helps make a visual proof of concept on arbitrary data, either synthetic or real. The code below does not yet do much, but once it’s finalised the code will be published properly.

Octave functions of relevance (still experimental and the GPLv3-licensed BioSig project may be required, as noted here):

function [new_image] = getdicomimages(directory, threshold)
% GETDICOMIMAGES: Open and test DICOM images of the heart
%
%    Licensing
%
%      GNU GPL version 3
%
%    GENERAL
%
%      Script - experimental, for handling MRI images with biosig
%
%    INPUT/S
%
%      -Omitted
%
%    OUTPUT/S
%
%      -Omitted
%
%    PENDING WORK
%
%      -
%
%    KNOWN BUG/S
%
%      -None.
%
%    COMMENT/S
%
%      -
%
%    RELATED FUNCTION/S
%
%      
%
%    ABOUT
%
%      -Created:     November 2010
%      -Last update: November 5th, 2010
%      -Revision:    0.0.1
%      -Author:      Dr. R. S. Schestowitz
% ==============================================================

%// addpath ("/home/roy/Main/IT/Programs/get_dicom_heart")
%// addpath ("/home/roy/Main/IT/Programs/biosig4octmat-2.50/biosig/t200_FileAccess")
%// img = opendicom("~/Main/IT/Programs/MR-Cardiac/i148713.MRDC.1");  // deprecated function
%// sopen("~/Main/IT/Programs/MR-Cardiac/Full-CD/DICOM/ST000000/SE000019/MR000000");
%// img(i)=imread('sombrero.jpg')
%// plot(img.T0)   //metadata attributes
%// ~/Main/IT/Programs/dicom2/dicom2 ~/Main/IT/Programs/MR-Cardiac/Full-CD/DICOM/ST000000/SE000014/*

cd(directory);
files_list=ls('*.tga');
n_images=size(files_list,1);
% BIOSIG open function, wrapper to DICOM handler:
% img = sopen("~/Main/IT/Programs/MR-Cardiac/i148713.MRDC.1");
for i=1:1:n_images,
	% files_list
	%// img_metadata(i) = sopen(files_list(i,:))
	%a=imread(files_list(i,:));
	%imshow(a)
	figure;
	img(:,:,i)=imread(files_list(i,:));
	imshow(img(:,:,i))
end

For actual testing:

function [new_image] = loadtaggedsequence(directory, threshold, draw_grid)
% LOADTAGGEDSQUENCE: Open and test tagged images of the heart
%
%    Licensing
%
%      GNU GPL version 3
%
%    GENERAL
%
%      Script
%
%    INPUT/S
%
%      -Omitted
%
%    OUTPUT/S
%
%      -directory: location of images,
%	 -threshold: colour of tags/markers based on conrtast of image
%      -draw_grid: boolean for display of manual grid/points (drawing or not drawing)
%
%    PENDING WORK
%
%      -
%
%    KNOWN BUG/S
%
%      -None.
%
%    COMMENT/S
%
%      -
%
%    RELATED FUNCTION/S
%
%      
%
%    ABOUT
%
%      -Created:     November 2010
%      -Last update: November 10th, 2010
%      -Revision:    0.0.1
%      -Author:      Dr. R. S. Schestowitz
% ==============================================================

% common commands

	%// loadtaggedsequence('~/Main/IT/Programs/MR-Cardiac/Zhen-thesis',1,1)
	%// addpath ("/home/roy/Main/IT/Programs/get_dicom_heart")
	%// addpath (genpath("/home/roy/Main/IT/Programs/biosig4octmat-2.50/biosig/"))
	%// path
	%// img = opendicom("~/Main/IT/Programs/MR-Cardiac/i148713.MRDC.1");  // deprecated function
	%// sopen("~/Main/IT/Programs/MR-Cardiac/Full-CD/DICOM/ST000000/SE000019/MR000000");
	%// img(i)=imread('sombrero.jpg')
	%// plot(img.T0)   //metadata attributes
	%// ~/Main/IT/Programs/dicom2/dicom2 ~/Main/IT/Programs/MR-Cardiac/Full-CD/DICOM/ST000000/SE000014/*

% experimental options

draw_circles=0;

% grid parameters

x_start = 10;
y_start = 10;
step_x = 30;
step_y = 20;
x_max = 110;
y_max = 80;


frame_size=8;
shuffle_radius=3;
 
% initialisation

first_image=1;
old_boundaries=[[0,0];[0,0]];

% get data

cd(directory);
files_list=ls('*.png');
	% files_list=ls('*.gif');
	% One image: files_list=('seq-tag-1.png');

% Opening [print filenames if necessary]

files_list
n_images=size(files_list,1);

% main loop

	% img = sopen("~/Main/IT/Programs/MR-Cardiac/i148713.MRDC.1");
for i=1:1:n_images,
		%// img_metadata(i) = sopen(files_list(i,:))
		%a=imread(files_list(i,:));
		%imshow(a)
	figure;
	%img(:,:,i)=imread(files_list(i,:));
	
	if (draw_grid==1)
		for j=x_start:step_x:x_max,
	   	  img(j,:,i)=0;
		end
		for j=y_start:step_y:y_max,
	   	  img(:,j,i)=0;
		end
	end
	if (draw_circles==1)
	     % [center,radius]=sphfit([10,11,0],[5,5,0]);	
        end
	boundaries=[[x_start,y_start];[x_start,y_start+step_y]; ...
	[x_start,y_start+step_y*2]; [x_start,y_start]; ...
	[x_start+step_x,y_start];[x_start+step_x*2,y_start]; ...
	[x_start+step_x*2,y_start+step_y];[x_start+step_x*2,y_start+step_y*2]; ...
	[x_start+step_x,y_start+step_y];[x_start+step_x,y_start+step_y*2]];
	%%	boundaries_x=[;    %% deprecated - no more xy separation
        	for l=1:1:size(boundaries,1)
        	%%	for k=1:1:size(boundaries,2)
        		% debugging
        		%	boundaries_y(k,1);
        		%	boundaries_y(2,l);
                     	% img(1,3,1)=1;
                     	% boundaries_y(k,1)
                     	% boundaries_y(2,l)
                     	%i
                     	%k
                     	%l
                     	% size(boundaries_y,1)
                     		if (first_image==0)
                     			%old_boundaries
                     	    		new_dot_x=old_boundaries(l,1)+1;
                     	    		new_dot_y=old_boundaries(l,2)+1;
                     	    			%old_boundaries(l,1)-8
                     	    		search_frame_current=img(((new_dot_x-frame_size):(new_dot_x+frame_size)),((new_dot_y-frame_size):(new_dot_y+frame_size))) 
                     	    		search_frame_previous=img(((old_boundaries(l,1)-frame_size):(old_boundaries(l,1)+frame_size)),((old_boundaries(l,2)-frame_size):(old_boundaries(l,2)+frame_size))) 
                     	    		shuffle_diff=shuffle_transform(search_frame_current, search_frame_previous, shuffle_radius)
                     	    		
                     	    		img(new_dot_x,new_dot_y,i)=threshold/255;
                     	    		old_boundaries(l,1)=new_dot_x;
                     	    		old_boundaries(l,2)=new_dot_y;
                     		else
                     	    		img(boundaries(l,1),boundaries(l,2),i)=(255-threshold)/255;
                     		end
                    %% 	end	
        	end
        	%for l=1:1:size(boundaries_x,1)
        		%for k=1:1:size(boundaries_x,2)
                     	   %img(boundaries_x(k,2),boundaries_x(l,1),i)=250;
                	%end
        	%end
	imshow(img(:,:,i))
		% plot3(img(:,:,i))
        if (first_image==1)
	   old_boundaries=boundaries;
	end
        first_image=0;
end

Peripheral transform:

function [value,shuffle] = shuffle_transform(image,im2,radius)
% SHUFFLE_TRANSFORM: A shuffle transform based on an implementation by Steve Marsland
%
%    Licensing
%
%      -
%
%    GENERAL
%
%      Function
%
%    INPUT/S
%
%      -image: first image
%      -im2: second image
%      -radius: the shuffle radius
%
%    OUTPUT/S
%
%      -Omitted
%
%    PENDING WORK
%
%      -
%
%    KNOWN BUG/S
%
%      -None.
%
%    COMMENT/S
%
%      -
%
%    RELATED FUNCTION/S
%
%      
%
%    ABOUT
%
%      -Created:     November 2010
%      -Last update: November 21st, 2010
%      -Revision:    0.0.1
%      -Author:      Dr. R. S. Schestowitz
% ==============================================================

shuffle = zeros(size(image));

for i=radius+1:size(image,1)-radius
    for j=radius+1:size(image,2)-radius
        template = im2(i-radius:i+radius,j-radius:j+radius);
        diff = abs(template-image(i,j));        
        shuffle(i,j) = min(diff(:));
    end
end
value = mean(shuffle(:))

Technical Notes About Comments

Comments may include corrections, additions, citations, expressions of consent or even disagreements. They should preferably remain on topic.

Moderation: All genuine comments will be added. If your comment does not appear immediately (a rarity), it awaits moderation as it contained a sensitive word or a URI.

Trackbacks: The URI to TrackBack this entry is:

https://schestowitz.com/Weblog/archives/2010/11/21/more-mri-tagging-papers/trackback/

Syndication: RSS feed for comments on this post RSS 2

    See also: What are feeds?, Local Feeds

Comments format: Line and paragraph breaks are automatic, E-mail address never displayed, HTML allowed: <a href="" title=""> <abbr title=""> <acronym title=""> <b> <blockquote cite=""> <cite> <code> <del datetime=""> <em> <i> <q cite=""> <strike> <strong>

Back to top

Retrieval statistics: 21 queries taking a total of 0.097 seconds • Please report low bandwidth using the feedback form
Original styles created by Ian Main (all acknowledgements) • PHP scripts and styles later modified by Roy Schestowitz • Help yourself to a GPL'd copy
|— Proudly powered by W o r d P r e s s — based on a heavily-hacked version 1.2.1 (Mingus) installation —|