phase_lines

status =phase_lines(freqs,obs,array,lines,phasefname);

freqs = frequencies to plot, e.g. [3 30] limits the plot at 3 and 30 Hz  
obs = observations in the power file to plot (defaults to all)
array = either an array name,e.g., 'all','hemisphere','quadrant','frontback',
    or 'oned' or a list of channels
 Restrictions nobs <= 3. If nobs == 1 and multiple frequencies 
 are listed then freqs < = 3. If nobs > 1 obs go on same plot one for
 each freq.  if nobs == 1 all the freqs go on the same plot
lines = 'black' or 'color'
phasefname = phasefilename (can be gui by skipping this param)

note also draws in solid black phase due to uncorrelated sources
in the four spheres model.  If 'laplacian' is a reference the curve
for laplacians is drawn as well.

0001 function status =phase_lines(freqs,obs,array,lines,phase_type,phasefname);
0002 %status =phase_lines(freqs,obs,array,lines,phasefname);
0003 %
0004 %freqs = frequencies to plot, e.g. [3 30] limits the plot at 3 and 30 Hz
0005 %obs = observations in the power file to plot (defaults to all)
0006 %array = either an array name,e.g., 'all','hemisphere','quadrant','frontback',
0007 %    or 'oned' or a list of channels
0008 % Restrictions nobs <= 3. If nobs == 1 and multiple frequencies
0009 % are listed then freqs < = 3. If nobs > 1 obs go on same plot one for
0010 % each freq.  if nobs == 1 all the freqs go on the same plot
0011 %lines = 'black' or 'color'
0012 %phasefname = phasefilename (can be gui by skipping this param)
0013 %
0014 %note also draws in solid black phase due to uncorrelated sources
0015 %in the four spheres model.  If 'laplacian' is a reference the curve
0016 %for laplacians is drawn as well.
0017 if nargin < 6
0018     [fid] = get_fid('rb');
0019 else
0020     fid = fopen(phasefname,'rb');
0021 end;
0022 version = fread(fid,1,'int16');
0023 if version ~= -2
0024     error('this is not a phase file');
0025 end;
0026  [nfiles,obs_labels,Epoch,Window_Length,NEpoch,Nbad_chan,bad_chan,ref_flag,reference,NChan,NFreq] = rd_anal_hdr(fid);
0027 if nargin == 0
0028     freqs = [0:40];
0029     obs = [1:nfiles];
0030     array = 'all';
0031     lines = 'black';
0032     phase_type = 'angles';
0033 end;
0034 if nargin == 1
0035     obs = [1:nfiles];
0036     array = 'all';
0037     lines = 'black';
0038     phase_type = 'angles';
0039 end;        
0040 if nargin == 2
0041     array = 'all';
0042     lines = 'black';
0043     phase_type = 'angles';
0044 end;    
0045 if nargin == 3
0046     lines = 'black';
0047     phase_type = 'angles';
0048 end;
0049 if nargin == 4
0050     phase_type = 'angles';
0051 end;
0052 if isempty(lines)
0053     lines = 'black';
0054 end
0055 if isempty(array)
0056     array = 'all';
0057 end
0058 if isempty(obs)
0059     obs = [1:nfiles];
0060 end;
0061 if isempty(freqs)
0062     freqs = [0:40];
0063 end;
0064 if isempty(phase_type)
0065     phase_type = 'angles';
0066 end;
0067 if size(obs,2) > 3
0068     fclose('all')
0069     error('too many observations > 3 to be plotted in a single plot')
0070 end;
0071 if size(obs,2) == 1
0072     if size(freqs,2) > 3 
0073         error('too many frequencies for a single plot')
0074     end;
0075 end;
0076 ch_pair_indices;
0077 if strcmp(lines,'black')
0078 line_type = ['k- ';'k--';'k-.'];
0079 scat_type = ['kx';'ko';'k.'];
0080 else    
0081 line_type = ['b-';'r-';'k-'];
0082 scat_type = ['kx';'ko';'k.'];
0083 end;
0084 
0085 iline = 1;
0086 [xelec,yelec,zelec] = electrodes(129);
0087 twod_dist = twod_pos(xelec,yelec,zelec,9.2);
0088 phase_array = real(get_eeg_coherence(fid,freqs,obs,nfiles,NChan,NFreq,Epoch));
0089 if strcmp(phase_type,'cos')
0090     phase_array = cos(phase_array);
0091     phasemax = 1;
0092 else
0093     phasemax = pi;
0094 end;
0095 if isstr(array)
0096     [arrayg,gnames] = arrays(array,[]);
0097 else
0098     arrayg = array;
0099 end
0100 if size(obs,2) > 1
0101     for g = 1:size(arrayg,1)
0102             for i= 1:size(freqs,2)
0103                 figure
0104                 hold on
0105                 
0106                 iline = 1;
0107                 
0108                 for j = 1:size(obs,2)
0109                 
0110                     if isstr(array)
0111                         [arrayh,hnames] = arrays(array,bad_chan(obs(j),find(bad_chan(obs(j),:))));
0112                     else
0113                         arrayh = array;
0114                     end
0115                 pairels = arrayh(g,find(arrayh(g,:)));
0116                 npairs = size(pairels,2)*(size(pairels,2)+1)/2;
0117                 pairs = zeros(1,npairs);
0118                 ipairs = 1;
0119                 for ip = 1:size(pairels,2)
0120                     for jp=ip:size(pairels,2)
0121                         pairs(ipairs) = chpair(pairels(ip),pairels(jp));
0122                         ipairs = ipairs + 1;
0123                     end;
0124                 end;
0125                 pairnames = hnames(g,:);
0126 %                hold on, plot(real(twod_dist(pairs)),phase_array((j-1)*size(freqs,2)+i,pairs),scat_type(iline,:))
0127                 p_coh = polyfit(twod_dist(pairs),phase_array((j-1)*size(freqs,2)+i,pairs),4);
0128                 coh_fit = real(polyval(p_coh,[0:1:fix(max(twod_dist(pairs)))]));
0129                 hold on, plot([0:1:max(twod_dist(pairs))],coh_fit,line_type(iline,:)), axis([0 max(twod_dist(pairs))+1 -phasemax*1.1 phasemax*1.1])
0130                 title(['Phase versus Distance:  @' int2str(freqs(i)) ' Hz ' hnames(g,:)]);
0131                 hold on, plot([max(twod_dist(pairs))-3  max(twod_dist(pairs))-2 ],[1-0.1*j 1-0.1*j],line_type(iline,:))
0132                 text([max(twod_dist(pairs))-2.5],[1-0.1*j],[deblank(setstr(obs_labels(obs(j),:)))])
0133                 xlabel('Interelectrode Distance (cm)')
0134                 ylabel('Phase')
0135                  iline = iline+ 1;
0136             end
0137         end
0138     end
0139 else
0140     if isstr(array)
0141         [arrayh,hnames] = arrays(array,bad_chan(obs,find(bad_chan(obs,:))));
0142     else
0143         arrayh = array;
0144     end
0145     for g = 1:size(arrayh,1)
0146         iline = 1;
0147         figure
0148         hold on
0149         
0150         for j = 1:size(freqs,2)
0151             pairels = arrayh(g,find(arrayh(g,:)));
0152             
0153             npairs = size(pairels,2)*(size(pairels,2)+1)/2;
0154             pairs = zeros(1,npairs);
0155             ipairs = 1;
0156             for ip = 1:size(pairels,2)
0157                 for jp=ip:size(pairels,2)
0158                     pairs(ipairs) = chpair(pairels(ip),pairels(jp));
0159                     ipairs = ipairs + 1;
0160                 end;
0161             end;
0162             chpair(pairels(1),pairels(2))
0163             pairnames = hnames(g,:);
0164             
0165 %            hold on, plot(real(twod_dist(pairs)),phase_array(j,pairs),scat_type(iline,:))
0166             p_coh = polyfit(real(twod_dist(pairs)),phase_array(j,pairs),4);
0167             coh_fit = real(polyval(p_coh,[0:1:fix(max(twod_dist(pairs)))]));
0168             hold on, plot([0:1:max(twod_dist(pairs))],coh_fit,line_type(iline,:)), axis([0 max(twod_dist(pairs))+1 -phasemax*1.1 phasemax*1.1])
0169             title(['Phase versus Distance: ' deblank(setstr(obs_labels(obs(1),:))) ' ' hnames(g,:)]);
0170             hold on, plot([max(twod_dist(pairs))-3 max(twod_dist(pairs))-2],[1-0.1*j 1-0.1*j],linetype(iline,:))
0171             text([max(twod_dist(pairs))-2.5],[1-0.1*j],[int2str(freqs(j)) ' Hz'])
0172             xlabel('Interelectrode Distance (cm)')
0173             ylabel('Phase')
0174             iline = iline+ 1;
0175         end
0176     end
0177 end
0178 status = 1;
0179 
0180 
0181 
0182 
0183 
0184 
0185 
0186 
0187 
0188 
0189 
0190 
0191 
0192 
0193 
0194 
0195 
0196 
0197 
0198 
0199 
0200 
0201 
0202 
0203 
0204 
0205 
0206 
0207