makeKymographC

 makeKymographC : creates a kymograph for given cell data file..
 A kymograph shows the fluorescence of the cell along the long axis
 of the cell, with time.

 INPUT :
       data : cell data file
       disp_flag : 1 to display image, 0 to not display image
       CONST : segmentation parameters
       which_channel : binarry array of fluorescence channels to be plotted eg. [1,1,1]

 OUTPUT :
       Kymo: Kymo has images at .r .g and .b fields. The (autogained)
       combination produces the kymgraph.
       ll1: size of the y axis in pixels.
       f1mm: is a two value array with the max and min value of channel 1
       f2mm: is a two value array with the max and min value of channel 2

 Copyright (C) 2016 Wiggins Lab 
 Written by Paul Wiggins, Stella Stylianidou.
 University of Washington, 2016
 This file is part of SuperSegger.
 
 SuperSegger is free software: you can redistribute it and/or modify
 it under the terms of the GNU General Public License as published by
 the Free Software Foundation, either version 3 of the License, or
 (at your option) any later version.
 
 SuperSegger is distributed in the hope that it will be useful,
 but WITHOUT ANY WARRANTY; without even the implied warranty of
 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 GNU General Public License for more details.
 
 You should have received a copy of the GNU General Public License
 along with SuperSegger.  If not, see <http://www.gnu.org/licenses/>.

0001 function [Kymo,ll1,f1mm,f2mm] = makeKymographC( data, disp_flag, CONST, which_channel )
0002 % makeKymographC : creates a kymograph for given cell data file..
0003 % A kymograph shows the fluorescence of the cell along the long axis
0004 % of the cell, with time.
0005 %
0006 % INPUT :
0007 %       data : cell data file
0008 %       disp_flag : 1 to display image, 0 to not display image
0009 %       CONST : segmentation parameters
0010 %       which_channel : binarry array of fluorescence channels to be plotted eg. [1,1,1]
0011 %
0012 % OUTPUT :
0013 %       Kymo: Kymo has images at .r .g and .b fields. The (autogained)
0014 %       combination produces the kymgraph.
0015 %       ll1: size of the y axis in pixels.
0016 %       f1mm: is a two value array with the max and min value of channel 1
0017 %       f2mm: is a two value array with the max and min value of channel 2
0018 %
0019 % Copyright (C) 2016 Wiggins Lab
0020 % Written by Paul Wiggins, Stella Stylianidou.
0021 % University of Washington, 2016
0022 % This file is part of SuperSegger.
0023 %
0024 % SuperSegger is free software: you can redistribute it and/or modify
0025 % it under the terms of the GNU General Public License as published by
0026 % the Free Software Foundation, either version 3 of the License, or
0027 % (at your option) any later version.
0028 %
0029 % SuperSegger is distributed in the hope that it will be useful,
0030 % but WITHOUT ANY WARRANTY; without even the implied warranty of
0031 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0032 % GNU General Public License for more details.
0033 %
0034 % You should have received a copy of the GNU General Public License
0035 % along with SuperSegger.  If not, see <http://www.gnu.org/licenses/>.
0036 
0037 white_bg = 0; % set to 1 for white background, outlined kymo
0038 Kymo = [];
0039 ll1=[];
0040 f1mm=[0,1];
0041 f2mm=[0,1];
0042 pixelsize = CONST.getLocusTracks.PixelSize;
0043 
0044 if isempty(pixelsize)
0045     pixelsize = 0.1;
0046 end
0047 
0048 
0049 if ~isfield(CONST.view, 'falseColorFlag' )
0050     CONST.view.falseColorFlag = false;
0051 end
0052 
0053 if ~exist( 'which_channel', 'var' ) || isempty(which_channel)
0054     which_channel = [1,1,1];
0055 end
0056 
0057 if ~isfield(CONST.view, 'filtered' )
0058     CONST.view.filtered = true;
0059 end
0060 
0061 filt_channel =  CONST.view.filtered ;
0062 
0063 
0064 
0065 persistent colormap_;
0066 if isempty( colormap_ )
0067     colormap_ = colormap( 'jet' );
0068     %colormap_ = colormap( 'hsv' );
0069 end
0070 
0071 if nargin < 2
0072     disp_flag = 0;
0073 end
0074 
0075 num_im = numel(data.CellA);
0076 
0077 ss = [0,0];
0078 ll = [0,0];
0079 
0080 % ss is max size of phase image
0081 % ll is max length/width of the cell
0082 for ii = 1:num_im
0083     ss_tmp = size(data.CellA{ii}.phase);
0084     
0085     ss(1) = max([ss(1),ss_tmp(1)]);
0086     ss(2) = max([ss(2),ss_tmp(2)]);
0087     
0088     ll_tmp = data.CellA{ii}.length;
0089     
0090     ll(1) = max([ll(1),ll_tmp(1)]);
0091     ll(2) = max([ll(2),ll_tmp(2)]);
0092 end
0093 
0094 ss = ss + [5,5];
0095 
0096 [XX,YY] = meshgrid( 1:ss(2), 1:ss(1) );
0097 ll0 = ll;
0098 ll = ceil(ll/2);
0099 
0100 ll1 = [-ll(1):ll(1)];
0101 ll2 = [-ll(2):ll(2)];
0102 
0103 [LL1,LL2] = meshgrid( ll1,ll2 );
0104 
0105 nn = numel( ll1 );
0106 
0107 kymoR = zeros(nn,num_im);
0108 kymoG = zeros(nn,num_im);
0109 kymoB = zeros(nn,num_im);
0110 
0111 
0112 for ii = 1:num_im
0113     
0114     mask  = data.CellA{ii}.mask;
0115     
0116     if isfield( data.CellA{ii}, 'fluor1') && which_channel(1)
0117         if filt_channel && isfield( data.CellA{ii},'fluor1_filtered')
0118             fluor1 =data.CellA{ii}.fluor1_filtered;
0119         else
0120             fluor1  = data.CellA{ii}.fluor1;
0121             if isfield( data.CellA{ii}, 'fl1' ) && isfield( data.CellA{ii}.fl1, 'bg' )
0122                 fluor1 = fluor1 - data.CellA{ii}.fl1.bg;
0123                 fluor1(fluor1<0) = 0;
0124             else
0125                 fluor1 = fluor1 - mean( fluor1(mask));
0126                 fluor1(fluor1<0) = 0;
0127             end
0128         end
0129     else
0130         fluor1 = 0*mask;
0131     end
0132     
0133     if isfield( data.CellA{ii}, 'fluor2') && which_channel(2)
0134         if filt_channel && isfield( data.CellA{ii}, 'fluor2_filtered' )
0135             fluor2 = data.CellA{ii}.fluor2_filtered;
0136         else
0137             fluor2 = data.CellA{ii}.fluor2;
0138             if isfield( data.CellA{ii}, 'fl12' ) && ...
0139                     isfield( data.CellA{ii}.fl2, 'bg' )
0140                 fluor2 = fluor2 - data.CellA{ii}.fl2.bg;
0141                 fluor2(fluor1<0) = 0;
0142             else
0143                 fluor2 = fluor2 - mean( fluor2(mask));
0144                 fluor2(fluor2<0) = 0;
0145             end
0146         end
0147     else
0148         fluor2 = 0*data.CellA{ii}.mask;
0149     end
0150     
0151     mask  = data.CellA{ii}.mask;
0152     
0153     % Make all the images the same sizes
0154     [rChan,~] = (fixIm(double((fluor2)).*double(mask),ss));
0155     [gChan,~] = (fixIm(double((fluor1)).*double(mask),ss));
0156     [bChan,roffset] = fixIm(mask,ss);
0157     
0158     ro = data.CellA{ii}.r_offset;
0159     r = data.CellA{ii}.r;
0160     
0161     e1 = data.CellA{ii}.coord.e1;
0162     e2 = data.CellA{ii}.coord.e2;
0163     
0164     LL1x =  LL1*e1(1)+LL2*e2(1)+r(1)-ro(1)+1+roffset(1);
0165     LL2y =  LL1*e1(2)+LL2*e2(2)+r(2)-ro(2)+1+roffset(2);
0166     
0167     rChanp = (interp2(XX,YY,double(rChan),LL1x,LL2y));
0168     gChanp = (interp2(XX,YY,double(gChan),LL1x,LL2y));
0169     bChanp = (interp2(XX,YY,double(bChan),LL1x,LL2y));
0170     
0171     rChanps = sum( double(rChanp) );
0172     gChanps = sum( double(gChanp) );
0173     bChanps = sum( double(bChanp) );
0174     
0175     kymoR(:,ii) = rChanps';
0176     kymoG(:,ii) = gChanps';
0177     kymoB(:,ii) = bChanps';
0178 end
0179 
0180 Kymo = [];
0181 
0182 if ~isfield(data.CellA{1}, 'pole');
0183     data.CellA{1}.pole.op_ori = 1;
0184 end
0185 
0186 
0187 if data.CellA{1}.pole.op_ori < 0 % flip the kymograph
0188     Kymo.g = kymoG(end:-1:1,:);
0189     Kymo.b = kymoB(end:-1:1,:);
0190     Kymo.b(isnan(Kymo.b)) = 0;
0191     Kymo.r = kymoR(end:-1:1,:);
0192 else
0193     Kymo.g = kymoG;
0194     Kymo.b = kymoB;
0195     Kymo.b(isnan(Kymo.b)) = 0;
0196     Kymo.r = kymoR;
0197 end
0198 
0199 f1mm(1) = min(Kymo.g(logical(Kymo.b)));
0200 f1mm(2) = max(Kymo.g(logical(Kymo.b)));
0201 f2mm(1) = min(Kymo.r(logical(Kymo.b)));
0202 f2mm(2) = max(Kymo.r(logical(Kymo.b)));
0203 
0204 
0205     if CONST.view.falseColorFlag
0206         % false color figure
0207         backer3 = double(cat(3, Kymo.b, Kymo.b, Kymo.b)>1);
0208         im = doColorMap( ag(Kymo.g,f1mm(1), f1mm(2)), colormap_ );
0209         imm =  im.*backer3+.6*(1-backer3);
0210     elseif white_bg     
0211         % figure with outline and white background
0212         sq = [1 1 1 ; 1 1 1 ; 1 1 1];
0213         backer = (ag(~Kymo.b));
0214         outline = imdilate(backer,sq) - backer;
0215         imm = cat(3, (ag(Kymo.r))+backer+0.2*outline, ...
0216             (ag(Kymo.g))+backer+0.2*outline,...
0217             backer+0.6*outline);        
0218     else
0219         % figure without outline and gray background
0220         backer = (ag(Kymo.b));
0221         backer = 0.3*(max(backer(:))-backer);
0222         imm = cat(3, (ag(Kymo.r))+backer, ...
0223             (ag(Kymo.g))+backer,...
0224             backer);       
0225     end
0226     
0227     if disp_flag
0228         figure (2);
0229         clf;
0230         imagesc(1:num_im,pixelsize*ll1, imm);
0231     end
0232 
0233 end
0234 
0235 function [imFix,roffset] = fixIm(im, ss)
0236 
0237 ssOld = size(im);
0238 imFix = zeros(ss);
0239 
0240 offset = floor((ss-ssOld)/2)-[1,1];
0241 if offset(1)<0
0242     offset(1) = offset(1) + 1;
0243 end
0244 if offset(2)<0
0245     offset(2) = offset(2) + 1;
0246 end
0247 
0248 try
0249     imFix(offset(1)+(1:ssOld(1)),offset(2)+(1:ssOld(2))) = im;
0250 catch ME
0251     printError(ME);
0252 end
0253 roffset = offset(2:-1:1);
0254 end
0255 
0256