Merge pull request #21 from openep/develop

Pull request for first release

README.md

@@ -0,0 +1,3 @@
+# OpenEP
+
+OpenEP supports the import and analysis of electroanatomic mapping data in Matlab

batchImport.m

@@ -0,0 +1,125 @@
+% ----------------------------------------------------------------------- %
+% OPENEP/batchImport is a template script to automate the import of 
+% datasets into OpenEP data format. This batchImport script requires a
+% number of assumptions to be true in order to function.
+%
+%   1. Cases should be named using a study number (e.g. 001) and compressed
+%   into folders named that study number.
+%   2. An Excell information spreadsheet should be provided which 
+%   identifies the map to import and should contain two columns where the 
+%   first contains the study number and the second the number of points per 
+%   map.
+%
+% As currently configured, the script assumes that datasets are stored on
+% an external hard drive, but copies dataset zip files to the local machine
+% hard drive prior to unzipping, which is usually fastest. 
+%
+% After importing into OpenEP dataformat, intermediate files which can
+% number in the 1000s and occupy Gb of space, are removed.
+% ----------------------------------------------------------------------- %
+
+% ----------------------------------------------------------------------- %
+%                             Configuration
+% testmode copies and unzips files but does not process them
+testmode = false; 
+
+% Set up working directories and paths
+zip_dir = '/Volumes/Extreme SSD/openep/openep_carto_export_zips';
+working_dir = '/Users/Steven/Desktop/openep_working_dir';
+info_spreadsheet = '/Volumes/Extreme SSD/openep/openep_datasheet.xlsx';
+
+% Cell array of possible start of names of the study XML file
+startStrings = {'CARTO' 'PVI' 'Study' '1008' 'PERS' 'AF' 'LA' 'PAF' 'AIAT'};
+% ----------------------------------------------------------------------- %
+
+% Load the study dataset Excel file
+info = xlsread(info_spreadsheet, 'A:B');
+
+% Get a list of files in zip_dir
+disp('Getting list of filenames')
+allFiles = nameFiles(zip_dir, 'showhiddenfiles', false);
+
+% Loop through all the zip files
+for i = 1:numel(allFiles)
+    disp(['Processing case: ' allFiles{i}]);
+
+    % Copy the zip file to the hard drive
+    sourceFile = [zip_dir filesep() allFiles{i}];
+    destinationFile = [working_dir filesep() allFiles{i}];
+    status = copyfile(sourceFile, destinationFile);
+
+    % Get the study number
+    studyNumber = str2double(allFiles{i}(1:3));
+
+    % Unzip the zip file to a folder named studyNumber in working_dir
+    caseDir = ([working_dir filesep() num2str(studyNumber)]);
+    unzip(destinationFile, caseDir);
+
+    % Locate the carto export xml file
+    allXmlFiles = nameFiles(caseDir, 'extension', 'xml', 'showhiddenfiles', false);
+    iF = NaN(size(startStrings));
+    for j = 1:numel(startStrings)
+       temp = find(strstartcmpi(startStrings{j}, allXmlFiles)); 
+       if ~isempty(temp)
+           iF(j) = temp;
+       end
+    end
+
+    if nnz(~isnan(iF)) > 1
+        % Multiple possible files found so the user has to choose which xml file to load
+        iF(isnan(iF)) = [];
+        beep()
+        disp([num2str(numel(iF)) ' possible XML files found:'])
+        for j = 1:numel(iF)
+           disp([ num2str(j) '    ' allXmlFiles{iF(j)}]);
+        end  
+        result = input('Enter choice: ');
+        studyXmlFile = [case_dir filesep() allXmlFiles{iF(result)}];
+    else % Automatically select the relevant xml file
+        iF(isnan(iF)) = [];
+        studyXmlFile = [caseDir filesep() allXmlFiles{iF}];
+    end
+
+    % Check that an XML file has been identified
+    if ~isfile(studyXmlFile)
+       error(['RUN_PEPR_EXPERIMENT: Specified XML file - ' studyXmlFile ' - does not exist']); 
+    end
+
+    % Read the number of points from the spreadsheet
+    numpts = info(info(:,1)==studyNumber,2);
+
+    % Specify the output directory and filename
+    outputFile = ([working_dir filesep() [num2str(studyNumber) '.mat']]);
+
+    % Start a timer
+    tic();
+
+    % Run importcarto_mem from the command line
+    if ~testmode
+        userdata = importcarto_mem(studyXmlFile ...
+            , 'maptoread', numpts ...
+            , 'refchannel', 'CS9-CS10' ...
+            , 'ecgchannel', 'V1' ...
+            , 'savefilename', outputFile ...
+            );
+    end
+
+    % Stop a timer
+    elapsedTime = toc();
+
+    % Get the size of caseDir and destinationFile
+    sizeUnZipped = du(['-sh ' caseDir]);
+    sizeZipped   = du(['-sh ' destinationFile]);
+
+    % Remove the carto files and zip file
+    deleteFolder(caseDir);
+    delete(destinationFile);
+
+    % Dump the time and filesizes to a text file
+    fid = fopen([working_dir filesep() [num2str(studyNumber) '.data']], 'wt');
+    fprintf(fid, [num2str(elapsedTime) '\n']);
+    fprintf(fid, [sizeUnZipped '\n']);
+    fprintf(fid, [sizeZipped '\n']);
+    fclose(fid);
+
+end

batchProcess.m

@@ -0,0 +1,52 @@
+% ----------------------------------------------------------------------- %
+% OPENEP/batchProcess is a template script which can be used to perform
+% multiple functions on OpenEP datasets
+% ----------------------------------------------------------------------- %
+
+% ----------------------------------------------------------------------- %
+%                             Configuration
+% Set up directories and load filenames
+working_dir = '/Users/Steven/Desktop/openep_working_dir';
+% ----------------------------------------------------------------------- %
+
+
+% Get a list of files
+disp('Getting list of filenames')
+allFiles = nameFiles(working_dir, 'showhiddenfiles', false, 'extension', 'mat');
+
+lat.carto = [];
+lat.openep = [];
+
+% Iterate through each file and perform some OpenEP functions
+for i = 1:numel(allFiles)
+   disp(['working on file ... ' allFiles{i}])
+   load([working_dir filesep() allFiles{i}]);
+
+   % Geometry
+   openep_num_pts(i,1) = getNumPts(userdata);
+   openep_chamber_area(i,1) = getArea(userdata);
+   openep_chamber_area_closed(i,1) = getArea(userdata, 'method', 'fill');
+   openep_chamber_volume(i,1) = getVolume(userdata);
+
+   % Local activation time
+   [~, earliestDirc(i,1:3)] = getEarliestActivationSite(userdata, 'method', 'ptbased');
+   [~, earliestPrct(i,1:3)] = getEarliestActivationSite(userdata, 'method', 'ptbasedprct');
+   tat.ptbased(i,1) = getTotalActivationTime(userdata, 'method', 'ptbased');
+   tat.ptbasedprct(i,1) = getTotalActivationTime(userdata, 'method', 'ptbasedprct');
+   tat.clinmap(i,1) = getTotalActivationTime(userdata, 'method', 'clinmap');
+   tat.clinmapprct(i,1) = getTotalActivationTime(userdata, 'method', 'clinmapprct');
+   tat.openepmap(i,1) = getTotalActivationTime(userdata, 'method', 'openepmap');
+   tat.openepmapprct(i,1) = getTotalActivationTime(userdata, 'method', 'openepmapprct');
+   lat_carto_temp = userdata.surface.act_bip(:,1);
+   lat_openep_tmp = generateInterpData(userdata, 'lat-map');
+   lat.carto = [lat.carto; lat_carto_temp];
+   lat.openep = [lat.openep; lat_openep_tmp];
+
+   % Voltage
+   meanVoltage.Carto(i,1) = getMeanVoltage(userdata, 'method', 'map');
+   meanVoltage.OpenEP(i,1) = getMeanVoltage(userdata, 'method', 'egm');
+   lva.Carto(i,1) = getLowVoltageArea(userdata, 'method', 'map');
+   lva.OpenEP(i,1) = getLowVoltageArea(userdata, 'method', 'egm');
+
+   % Conduction Velocity
+end

cartoMesh2VTK.m

@@ -0,0 +1,32 @@
+function tr = cartoMesh2VTK()
+% CARTOMESH2VTK Converts a Carto mesh file to VTK file
+%
+% Usage:
+%   tr = cartoMesh2VTK('openfile')
+% Where:
+%   tr,  - a TriRep object
+%
+% CARTOMESH2VTK Converts a Carto3 mesh to a VTK file and returns a TriRep
+% object
+%
+% Author: Steven Williams (2015) (Copyright)
+% SPDX-License-Identifier: Apache-2.0
+%
+% Modifications -
+%
+% Info on Code Testing:
+% ---------------------------------------------------------------
+% test code
+% ---------------------------------------------------------------
+%
+% ---------------------------------------------------------------
+% code
+% ---------------------------------------------------------------
+
+% output the mesh to a temporary VTK
+[fileName, pathName] = uigetfile('*.mesh', 'Select a Carto3 .mesh file');
+tr = read_meshfile([pathName filesep fileName]);
+fileName2 = [fileName(1:end-4) 'vtk'];
+writeTriRep2VTK(repack(tr), [], 'outputfile', [pathName filesep fileName2]);
+
+end

cemrg2carto.m

@@ -0,0 +1,58 @@
+function cemrg2carto(infile, outfile, vtktitle)
+% CEMRG2CARTO Converts a VTK for loading into Carto
+% Usage:
+%   op = cemrg2carto(infile, outfile)
+%
+% Where:
+%   infile - a VTK file created by CEMRG
+%   outfile - a VTK file formatted for Carto
+%
+% CEMRG2CARTO Converts a VTK file created by CEMRG into a VTK file that can
+% be loaded into Carto. The two main functions are to convert cell data to
+% point data and to add the correct header based on patient name and ID.
+%
+% Author: Steven Williams (2016)
+% Modifications -
+%
+% Info on Code Testing:
+% ---------------------------------------------------------------
+% 
+% ---------------------------------------------------------------
+%
+% ---------------------------------------------------------------
+% code
+% ---------------------------------------------------------------
+
+if strcmpi(infile, 'openfile')
+    [filename, pathname] = uigetfile('*.vtk', 'Select VTK file to open ...');
+    infile = [pathname filename];
+end
+if infile==0
+    disp('Operation cancelled');
+    return
+end
+
+h = waitbar(0, 'Load VTK file');
+hVtk = VTKReader(infile);
+
+waitbar(.2, h, 'Read VTK data');
+hVtk.readAllData();
+
+waitbar(.5, h, 'Get anatomy from VTK file');
+tr = hVtk.getTriRep();
+
+waitbar(.8, h, 'Convert cell data to point data');
+pointdata = trFaceToVertData(tr, hVtk.CellData{1});
+close(h);
+
+if exist('vtktitle')==0
+    answer = inputdlg({'firstname', 'surname', 'ID'});
+    vtktitle = ['PatientData ' answer{1} ' ' answer{2} ' ' answer{3}];
+end
+
+if strcmpi(outfile, 'openfile')
+    [filename, pathname] = uiputfile('*.vtk', 'Select VTK file to save ...');
+    outfile = [pathname filename];
+end
+
+writeTriRep2VTK(tr, pointdata, 'title', vtktitle, 'outputfile', outfile);

cvHistogram.m

@@ -0,0 +1,59 @@
+function cvHistogram( userdata, varargin )
+% CVHISTOGRAM Draws a conduction velocity histogram
+%
+% Usage:
+%   cvHistogram( userdata )
+% Where:
+%   userdata  - see importcarto_mem
+%
+% CVHISTOGRAM accepts the following parameter-value pairs
+%   'limits'    {[0 5]} | array
+%   'binwidth'  {0.1} | double
+%
+% CVHISTOGRAM displays a histogram of conduction velocities. Limits are set
+% to exclude non-physiological conduction velocities
+%
+% Author: Steven Williams (2020) (Copyright)
+% SPDX-License-Identifier: Apache-2.0
+%
+% Modifications -
+%
+% Info on Code Testing:
+% ---------------------------------------------------------------
+% cvHistogram( userdata )
+% ---------------------------------------------------------------
+%
+% ---------------------------------------------------------------
+% code
+% ---------------------------------------------------------------
+
+nStandardArgs = 1; % UPDATE VALUE
+limits = [0 5];
+binwidth = 0.1;
+if nargin > nStandardArgs
+    for i = 1:2:nargin-nStandardArgs
+        switch varargin{i}
+            case 'limits'
+                limits = varargin{i+1};
+            case 'binwidth'
+                binwidth = varargin{i+1};
+        end
+    end
+end
+% TODO: check format of input values for each parameter are correct
+
+cvdata = getConductionVelocity(userdata);
+cvdata(cvdata<limits(1)) = [];
+cvdata(cvdata>limits(2)) = [];
+
+histogram(cvdata ...
+    , 'normalization', 'probability' ...
+    , 'binwidth', binwidth ...
+    , 'facecolor', 'k' ...
+    , 'facealpha', 1 ...
+    );
+set(gcf, 'color', 'w');
+xlabel('Conduction Velocity (m/s)');
+ylabel('Frequency')
+set(gca, 'fontsize', 16)
+end