%Find input parameter ranges for central case %and Convert to natural units clear all ; N=5000; halfN=N/2; density=15000; % 1/micron**2 h=0.05 ; % microns %central case TfCenter=1207 %microseconds AmpCenter=873 %channels serial=(1:1568)'; %=====The normal data========================= fname='norm_avg_blip_results1.txt' fullname = ( strcat('c:\My Documents\BigData\bartol\april2000runs\May22format\',fname) ); [str1,str2,tf,ba,R,kp,D,peak,tr20,tr80,tpeak,tf90,tf33] = ... textread(fullname,... '%s%s%d%d%f%f%f%f%f%f%f%f%f',... 1568); %kp is in units of 1e7 1/Mole sec %D is in units of 1e-6 cm**2/sec %tf is in mSec alpha = (1+R)./ (tf*.001); beta = ba .* alpha; km = beta ./(2*R); %Tkm= 1e6 ./km ; %Talphabeta = 1e6 ./(alpha+beta); Tbind = 1.39 * h ./ (density * kp*1e7 * 1e15/6.02e23) *1e6; Tdiff = 0.8/(4*pi) * N ./(density*D*1e-6*1e8) *1e6 ; %symmetry = (tpeak - tr80)./(tf90 - tr80); peakNorm = peak; eff = peak ./ halfN ; risetime = tr80 - tr20; tr20norm=3*tr20./risetime; falltime = tf33 - tf90; plateautime = (tf90 - tr80)./sqrt(risetime.^2+falltime.^2); %=====The DFP data========================= fname='dfp_avg_blip_results.txt' fullname = ( strcat('c:\My Documents\BigData\bartol\april2000runs\May22format\',fname) ); [str1,str2,tf,ba,R,kp,D,peak,tr20,tr80,tpeak,tf90,tf33] = ... textread(fullname,... '%s%s%d%d%f%f%f%f%f%f%f%f%f',... 1568); peakDFP = peak; effDFP = peak ./ halfN ; risetimeDFP = tr80 - tr20; falltimeDFP = tf33 - tf90; tr20DFP=3*tr20./risetimeDFP; plateautimeDFP = (tf90 - tr80)./sqrt(risetimeDFP.^2+falltimeDFP.^2); %=====The BTX data========================= fname='hbtx_avg_blip_results.txt' fullname = ( strcat('c:\My Documents\BigData\bartol\april2000runs\May22format\',fname) ); [str1,str2,tf,ba,R,kp,D,peak,tr20,tr80,tpeak,tf90,tf33] = ... textread(fullname,... '%s%s%d%d%f%f%f%f%f%f%f%f%f',... 1568); peakBTX=peak; effBTX = peak ./ halfN ; risetimeBTX = tr80 - tr20; falltimeBTX = tf33 - tf90; tr20BTX=3*tr20./risetimeBTX; plateautimeBTX = (tf90 - tr80)./sqrt(risetimeBTX.^2+falltimeBTX.^2); %=====Do conversions to real units=================== ampScale = peakNorm ./ AmpCenter; timeScale = falltime ./ TfCenter; Nreal = N ./ ampScale; alpha = alpha .* timeScale; beta = beta .* timeScale; Kp = kp .* timeScale; Km = km .* timeScale; Dreal = D .* timeScale./ampScale; peakNorm = peakNorm./ampScale; peakDFP = peakDFP./ampScale; peakBTX = peakBTX./ampScale; risetime = risetime ./ timeScale; risetimeDFP = risetimeDFP ./ timeScale; risetimeBTX = risetimeBTX ./ timeScale; falltime = falltime ./ timeScale; falltimeDFP = falltimeDFP ./ timeScale; falltimeBTX = falltimeBTX ./ timeScale; %=====Compute the acceptance ranges======== TargetRange = 0.2; % +/- log units for Target1=[2.6 ];%2.1 3.1] for Target2=[1.16 ];%1.70 0.70] for Target3=[0.69 ];%1.10 0.20] %=====The output data====================== data = [ ba, R, kp, D,Tbind, Tdiff... beta, alpha, Kp, Km, Dreal, Nreal,... risetime, .001*falltime, peakNorm,plateautime, tr20norm, eff, ... risetimeDFP, .001*falltimeDFP, peakDFP, plateautimeDFP, tr20DFP, effDFP, ... risetimeBTX, .001*falltimeBTX, peakBTX, plateautimeBTX, tr20BTX, effBTX, ... log(risetime ./ falltime), ... log(falltimeDFP ./ falltime), ... log(risetimeBTX ./ risetimeDFP) ... zeros(1568,1),... serial,... ]; %eliminate rows with eff<.05 data = data(find(data(:,15)>=0.05*AmpCenter),:); %get rows in Tr/Tf acceptance range data = data(find(-data(:,31)<(Target1+TargetRange)),:); data = data(find(-data(:,31)>(Target1-TargetRange)),:); %get rows in TfDFP/TfNorm acceptance range data = data(find(data(:,32)<(Target2+TargetRange)),:); data = data(find(data(:,32)>(Target2-TargetRange)),:); %get rows in TrBTX/TrDFP acceptance range data = data(find(data(:,33)<(Target3+TargetRange)),:); data = data(find(data(:,33)>(Target3-TargetRange)),:); %Compute sum(abs(log(target))) data(:,34) = sqrt( ... (data(:,31)+Target1).^2 ... + (data(:,32)-Target2).^2 ... +(data(:,33)-Target3).^2 ); %sort by first 4 columns sortedData=sortrows(data,[34]); %normalize each parameter (column) to the range 0-1.0 permute = [1 2 3 4]; % 2 1 3 4 % 1 4 2 3 for i=1:length(permute) pos(:,permute(i)) = .1*rand(size(data(:,i))) + ... (log(data(:,i))-log(min(data(:,i))))... /(log(max(data(:,i)))-log(min(data(:,i)))); end pos = pos'; npoints = length(data) ; mapres = 10 ; map = jet(mapres) ; %define a color map jet dep=34; %index of the dependent variable to plot s = (mapres)*... (data(:,dep)-min(data(:,dep)))/((max(data(:,dep))-min(data(:,dep))))+1 ; %s = 1:npoints; highbrk = npoints; lowbrk = 1; %set(gcf,'position',[300 0 600 800] ) %make a bigger figure %orient tall %fill the paper when we print cindex = fix(s(lowbrk)); plot(pos(:,lowbrk), 'color', 'white' ); %make the axes axisx = [1,4; 1,1; 2,2; 3,3; 4,4 ]'; axisy = [0,0; 0,1; 0,1; 0,1; 0,1 ]'; line(axisx, axisy, 'color', 'k','erasemode','none'); drawnow; axisx=[1,2 ; 2,3; 3,4 ]'; set(gca, 'drawmode','fast') for i = highbrk:-1:lowbrk if (s(i)<12) cindex = max(1,min(fix(s(i)),mapres)) ; axisy=[pos(1,i),pos(2,i); pos(2,i),pos(3,i);... pos(3,i),pos(4,i)]'; line(axisx,axisy,'erasemode','none','color',map(cindex,:),'linewidth',2); drawnow end end end end end %colorbar title ('Color is distance to target -- blue=close') set(gca,'xtick',[1 2 3 4]) set(gca,'xticklabel',['b/a|b/2Km|Kp|D']) ylabel('normalized parameter range')