Description of trans_coord

0001 function    [pick12,pick22,samp] = ...
0002             trans_coord_sbi(meg_dir, sysid, pick1, pick2, samp)
0003 % Input
0004 %
0005 % Output
0006 % samp,pick12,pick22
0007 %
0008 % Copyright (C) 2011, ATR All Rights Reserved.
0009 % License : New BSD License(see VBMEG_LICENSE.txt)
0010 
0011 % --- convert cordinate from DICOM-shift or Dewar to DICOM
0012 %
0013 
0014 global vbmeg_inst; 
0015 define = vbmeg_inst.const; 
0016 
0017 head_dir = '../head/';    % MEGディレクトリからheaderディレクトリへの相対path
0018 pos_dir  = '../pos/';    % MEGディレクトリからposディレクトリへの相対path
0019 mri_dir  = '../mri/';    % MEGディレクトリからmriディレクトリへの相対path
0020 
0021 %
0022 % --- read SBI header
0023 %
0024 fname = sprintf('%s%s%s%s',meg_dir,head_dir,sysid,define.HEADER_EXTENSION);
0025 fid   = fopen(fname ,'r','l');
0026 if fid==-1, disp('XXX Error : Cannot open head file !!\n'); return; end
0027 
0028 % file info.
0029 H_file=meg_hread('file',fid);
0030 
0031 %
0032 % Read Positioning info
0033 %
0034 
0035 % positioning coil
0036 H_pos = meg_hread('pos',fid,H_file.i_pos_begin);
0037 npos  = H_pos.i_ncoil;  % number of positioning coils
0038 
0039 H_pos_rept      = meg_hread('pos_rept',fid,H_file.i_pos_rept_begin);
0040 pos_lcz_begin = H_pos_rept.i_data_begin;
0041 
0042 % MRI
0043 H_mri  = meg_hread('mri',fid,H_file.i_mri_begin);
0044 mark   = H_mri.st_axes(1:npos,:);  % MRI marker position on DICOM-shift cord.
0045 mrisys = H_mri.c_sysid{1};  % MRI system ID
0046 
0047 % LCZ
0048 H_lcz  = meg_hread('lcz',fid,H_file.i_lcz_begin);
0049 sphr   = H_lcz.f_model_r;  % radias of sphere model
0050 sph0   = H_lcz.st_model;  % origin of sphere model on DICOM-shift cord.
0051 
0052 fclose(fid); 
0053 
0054 %
0055 % Read Vivid pos file
0056 %
0057 
0058 % pos. data
0059 fname = sprintf('%s%s%s%s',meg_dir,pos_dir,sysid,define.POS_EXTENSION);
0060 fid   = fopen(fname ,'r','l');
0061 
0062 if fid==-1, disp('XXX Error : Cannot open pos file !!\n'); return; end
0063 
0064 fseek(fid,pos_lcz_begin,'bof');
0065 pos=zeros(npos,3);
0066 
0067 for i=1:npos
0068     a=fread(fid,8,'float');
0069     b=reshape(a,1,8);
0070     pos(i,:)=b(1:3);  % position of positioning coils on DEWAR cord.
0071 end
0072 
0073 fclose(fid); clear a b;
0074 
0075 
0076 % DICOM-shift
0077 
0078 % Read mri-table file
0079 fname1 = sprintf('%s%s%s/%s%s',...
0080     meg_dir, mri_dir, mrisys, mrisys, define.DICOM1_EXTENSION);
0081 
0082 fid=fopen(fname1);
0083 
0084 if fid==-1, 
0085     disp('XXX Error : Cannot open mri-table file !!\n'); 
0086     disp(sprintf('mri-table file=%s\n', fname1));
0087     return; 
0088 end
0089 
0090 a=fscanf(fid,'%s'); 
0091 fclose(fid);
0092 
0093 sx=str2num(a(findstr(a,'image_x_vector=')+15:findstr(a,'image_y_vector=')-1));
0094 sy=str2num(a(findstr(a,'image_y_vector=')+15:findstr(a,'image_z_vector=')-1));
0095 sz=str2num(a(findstr(a,'image_z_vector=')+15:findstr(a,'pixel_spacing=')-1));
0096 s0=str2num(a(findstr(a,'z0_slice_origin=')+16:findstr(a,'magnification=')-1));
0097 ss=str2num(a(findstr(a,'pixel_spacing=')+14:findstr(a,'z0_slice_origin=')-1));
0098 
0099 sp1=findstr(a,'slice_pitch=')+12; 
0100 sp2=findstr(a(sp1:size(a,2)),'te=')-1; 
0101 sp=str2num(a(sp1:sp1+sp2-1));
0102 nsl=str2num(a(findstr(a,'number_of_slices=')+17:findstr(a,'number_of_echos=')-1));
0103 
0104 % Read mri-table file
0105 fid=fopen(fname1);
0106 
0107 b=1; 
0108 while b,
0109     if ~isempty(findstr(char(fgets(fid)),'slice_location_table')), b=0; end
0110 end
0111 
0112 sloc = zeros(nsl,1);
0113 
0114 for i=1:nsl, 
0115     sloc(i)=str2num(char(fgets(fid))); 
0116 end
0117 fclose(fid);
0118 
0119 scorner = sz*min(sloc)+s0;
0120 dcms0   = 128*sx*ss+127*sy*ss+fix(nsl/2)*sz*sp+scorner; 
0121 dcms0   = dcms0.*1e-3;
0122 
0123 % from DICOM-shift to DICOM
0124 
0125 % MRI marker position on DICOM cord.
0126 mark2 = mark + repmat(dcms0,npos,1);  
0127 
0128 if sphr<0, 
0129     sph02=sph0+dcms0';  % origin of sphere model on DICOM cord.
0130 else, 
0131     sph02=[0,0,0]; disp('XXX : No suport, because not using VIVID'); 
0132     return; 
0133 end
0134 
0135 % from Dewar to DICOM
0136 nch = size(pick1,1);    % # of sensor
0137 
0138 pick1(:,4)=mod(pick1(:,4),pi); 
0139 pick1(:,5:6)=mod(pick1(:,5:6),2*pi); % normalize 0<pi or 2pi
0140 pick2(:,4)=mod(pick2(:,4),pi); 
0141 pick2(:,5:6)=mod(pick2(:,5:6),2*pi); % normalize 0<pi or 2pi
0142 
0143 % angle -> normal vector
0144 pick1(:,4:6)=[sin(pick1(:,4)).*cos(pick1(:,5)) ...
0145               sin(pick1(:,4)).*sin(pick1(:,5)) cos(pick1(:,4))]; 
0146 % angle -> normal vector
0147 pick2(:,4:6)=[sin(pick2(:,4)).*cos(pick2(:,5)) ...
0148               sin(pick2(:,4)).*sin(pick2(:,5)) cos(pick2(:,4))]; 
0149 
0150 % minimum distance between a positioning coil and detection coils
0151 for i=1:npos  
0152     sfromc(i) = min(sqrt(sum( ...
0153                 ( pick1(:,1:3) - repmat(pos(i,:),size(pick1,1),1) ).^2,2)));
0154 end
0155 
0156 % select 3 coils
0157 sel=find(sfromc<max(sfromc)); 
0158 
0159 % axes of Head cord. on Dewar cord.
0160 ex=pos(sel(3),:)-pos(sel(1),:); 
0161 ex=ex/norm(ex);
0162 ez=cross(ex,pos(sel(2),:)-pos(sel(1),:)); 
0163 ez=ez/norm(ez);
0164 ey=cross(ez,ex); 
0165 mc=[ex; ey; ez]; 
0166 
0167 % origin of Head cord. on Dewar cord.
0168 oc=sum((pos(sel(2),:)-pos(sel(1),:)).*ex)*ex + pos(sel(1),:); 
0169 
0170 % axes of Head cord. on DICOM cord.
0171 ex=mark2(sel(3),:)-mark2(sel(1),:); 
0172 ex=ex/norm(ex);
0173 ez=cross(ex,mark2(sel(2),:)-mark2(sel(1),:)); 
0174 ez=ez/norm(ez);
0175 ey=cross(ez,ex); 
0176 mm=[ex; ey; ez]; 
0177 
0178 % origin of Head cord. on DICOM cord.
0179 om=sum((mark2(sel(2),:)-mark2(sel(1),:)).*ex)*ex + mark2(sel(1),:); 
0180 
0181 trans = mm'*mc;
0182 
0183 % position of detection coils on DICOM cord.
0184 pick12(:,1:3)=(trans*(pick1(:,1:3)-repmat(oc,nch,1))')'+repmat(om,nch,1);  
0185 % position of compensation coils on DICOM cord.
0186 pick22(:,1:3)=(trans*(pick2(:,1:3)-repmat(oc,nch,1))')'+repmat(om,nch,1);  
0187 
0188 % direction of detection coils on DICOM cord.
0189 pick12(:,4:6)=(trans*pick1(:,4:6)')'; 
0190 % direction of compensation coils on DICOM cord.
0191 pick22(:,4:6)=(trans*pick2(:,4:6)')'; 
0192 
0193 samp.sph0   = sph02;        % origin of sphere model on DICOM cord.
0194 samp.sphr   = sphr;         % radius of sphere model [m]
0195 samp.mrisys = mrisys;       % MRI system ID
trans_coord_sbi

PURPOSE

SYNOPSIS

DESCRIPTION

CROSS-REFERENCE INFORMATION

SOURCE CODE
