diff --git a/gamma_bomb_2/GammaBomb_2_UserGuide.pdf b/gamma_bomb_2/GammaBomb_2_UserGuide.pdf deleted file mode 100644 index d49ad9a..0000000 Binary files a/gamma_bomb_2/GammaBomb_2_UserGuide.pdf and /dev/null differ diff --git a/gamma_bomb_2/code/Repair_LST_HDR.m b/gamma_bomb_2/code/Repair_LST_HDR.m deleted file mode 100644 index be80a48..0000000 --- a/gamma_bomb_2/code/Repair_LST_HDR.m +++ /dev/null @@ -1,77 +0,0 @@ -function Repair_LST_HDR(lst_files) -%% Function To Repair Lst mode header files if missing -%% Written by Daniel Chonde -%% Req: cell of list file names -%% Output: *.lst.hdr files when missing -%% Written: 3/31/2011 Updated: 11/14/2011 -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%This function repairs missing .lst.hdr files which are necessary for -%gettting the timing information for reconstruction. -%% - - -%%Repair LST.HDR files -if size(dir('*.lst'),1)~=size(dir('*.lst.hdr'),1) - disp('It appears that some LST header files are missing we will attempt to approximate them'); - for m=1:size(lst_files,1) - lst_hdr_files(m,:)=strcat(lst_files(m),'.hdr'); - end - clear m - - for m=1:size(lst_hdr_files,1) - if exist(lst_hdr_files{m,:},'file')==0 - if exist(regexprep(strtrim(lst_hdr_files{m,:}),'.lst',''),'file') - disp([regexprep(strtrim(lst_hdr_files{m,:}),'.lst','') ' found...ask Larry to make .lst.hdr files again...']) - copyfile(regexprep(strtrim(lst_hdr_files{m,:}),'.lst',''),lst_hdr_files{m,:}) - else - disp([lst_hdr_files{m,:} ' is missing...approximating...']) - %figure out creation time of log file and subtract 58 seconds - no_ext_name=regexp(lst_hdr_files{m,:},'.lst.hdr','split'); - %no_ext_name=regexp(lst_hdr_files{m-1,:},'.lst.hdr','split'); - log_file=strcat(no_ext_name{1},'.lst'); - log=dir(log_file); - log_date=log.date; - log_creation_date_time=regexp(log_date,' ','split'); - log_h_m_s=regexp(cell2mat(log_creation_date_time(1,2)),':','split'); - log_file_creation_time=str2double(log_h_m_s{1,1})*3600+str2double(log_h_m_s{1,2})*60+str2double(log_h_m_s{1,3}); - %lst_file_creation_time=log_file_creation_time-58; - start=calc_lst_duration(log_file); - lst_file_creation_time=log_file_creation_time-start; - new_hdr_hh=num2str(floor(lst_file_creation_time/3600)); - if size(new_hdr_hh,2)==1 - new_hdr_hh=strcat('0',new_hdr_hh); - end - new_hdr_mm=num2str(floor(lst_file_creation_time/60-floor((lst_file_creation_time/3600))*60)); - if size(new_hdr_mm,2)==1 - new_hdr_mm=strcat('0',new_hdr_mm); - end - new_hdr_ss=num2str(lst_file_creation_time-str2num(new_hdr_mm)*60-str2num(new_hdr_hh)*3600); - if size(new_hdr_ss,2)==1 - new_hdr_ss=strcat('0',new_hdr_ss); - end - - if m==1 - copyfile(fullfile('D:\Users\Dan\Codes\matlab\externals','template.lst.hdr'),lst_hdr_files{m,:}); - else - copyfile(lst_hdr_files{m-1,:},lst_hdr_files{m,:}); - end - fileattrib(lst_hdr_files{m,:},'+w'); - %fid=fopen(lst_hdr_files{m-1,:}); - fid=fopen(lst_hdr_files{m,:}); - F = fread(fid, '*char')'; - fclose(fid); - C=regexp(F(1,:),'\(hh:mm:ss)\s:=','split'); - F((size(C{1},2)+15):(size(C{1},2)+16))=new_hdr_hh; - F((size(C{1},2)+18):(size(C{1},2)+19))=new_hdr_mm; - F((size(C{1},2)+21):(size(C{1},2)+22))=new_hdr_ss; - fid2=fopen(lst_hdr_files{m,:},'w'); - fwrite(fid2,F); - fclose(fid2); - fileattrib(lst_hdr_files{m,:},'-w'); - end - end - clear C F fid fid2 log log_creation_date_time log_date log_file log_file_creation_time log_h_m_s lst_file_creation_time new_hdr_hh new_hdr_mm new_hdr_ss no_ext_name - end - clear m lst_hdr_files -end -end \ No newline at end of file diff --git a/gamma_bomb_2/code/TAC_process_data.m b/gamma_bomb_2/code/TAC_process_data.m deleted file mode 100644 index 99ab3f7..0000000 --- a/gamma_bomb_2/code/TAC_process_data.m +++ /dev/null @@ -1,18 +0,0 @@ -function TAC_process_data(Sum_Output_name) -if ~exist(fullfile(pwd,'TAC_process_data.bat'),'file') - copyfile('D:\Users\Dan\Codes\bat\TAC_process_data.bat') -end -if ~exist(fullfile(pwd,'TAC_process_data_cleanup.bat'),'file') - copyfile('D:\Users\Dan\Codes\bat\TAC_process_data_cleanup.bat') -end - -dos(['hist_PET -l ' Sum_Output_name '.lst -s ' Sum_Output_name '.s -L ' Sum_Output_name '.lor -A c:\bin\lorlookup_p256.dat -C 19,128']); -copyfile(strcat(Sum_Output_name,'.s.hc'),strcat(Sum_Output_name,'_hp.s.hc'),'f'); -delete(strcat(Sum_Output_name,'_p.lor*'),strcat(Sum_Output_name,'_d.lor*'),strcat(Sum_Output_name,'.s*'),strcat(Sum_Output_name,'_*r*.s*'),strcat(Sum_Output_name,'.*map')) -process_data_command_line=['TAC_process_data ' Sum_Output_name]; -dos_ex(process_data_command_line); -process_data_cleanup_command_line=['TAC_process_data_cleanup ' Sum_Output_name]; -dos(process_data_cleanup_command_line); -convertshc(strcat(Sum_Output_name,'.s.hc'),strcat(Sum_Output_name,'_hp.s.hc')) -copyfile(strcat(Sum_Output_name,'.s.hc'),strcat(Sum_Output_name,'_head_curve.s.hc'),'f'); -end \ No newline at end of file diff --git a/gamma_bomb_2/code/TAC_sort_lst_mode.m b/gamma_bomb_2/code/TAC_sort_lst_mode.m deleted file mode 100644 index 0dc2045..0000000 --- a/gamma_bomb_2/code/TAC_sort_lst_mode.m +++ /dev/null @@ -1,122 +0,0 @@ -function [Scanner_Duration, Series_Duration, frame_start_time, frame_stop_time, acquisition_start, lst_files, lst_files_whole]=TAC_sort_lst_mode(refpath,calc_t_rez) - -if nargin<2, calc_t_rez=1000;end - -cd(refpath) -lst_files=dir('*.lst'); -lst_files_whole=dir('*.lst'); -%order Frames -for m=1:size(lst_files,1) - File{m,1}=lst_files(m).name; -end -clear m -for p=1:size(File,1) - order(p,:)=regexp(cell2mat(File(p)),'.lst','split'); - order1(p,:)=mat2cell(order{p,1},size(order{p,1},1),size(order{p,1},2)); - order2(p,:)=regexp(order1{p,1},'Frame','split'); - if size(order2,2)==2 - order3(p,:)=[str2double(order2{p,2}) p]; - else - order3(p,:)=[str2double(order2{p,1}) p]; - end -end -File_unsor=File; -clear File -sortorder=sortrows(order3,1); -File=File_unsor(sortorder(:,2)); -lst_files_whole=lst_files_whole(sortorder(:,2)); -clear p -clear File_unsor -lst_files=File; -clear File - -%% Determine end time from file time information -%%%%%This code is the old way of determining time information, from reading -%%%%%frame time information and hdr file. A faster implementation is to -%%%%%just read the first frame time information and then the time marks. -% for m=1:size(lst_files_whole,1) -% frame_date{m,1}=lst_files_whole(m).date; -% frame_end_date_time(m,:)=regexp(cell2mat(frame_date(m)),' ','split'); -% h_m_s(m,:)=regexp(cell2mat(frame_end_date_time(m,2)),':','split'); -% frame_stop_time(m,1)=str2double(h_m_s{m,1})*3600+str2double(h_m_s{m,2})*60+str2double(h_m_s{m,3}); -% end -% clear m -% -% %%Determine start time from lst.hdr information - Repair_LST_HDR(lst_files) -% for m=1:size(lst_files_whole,1) -% fid = fopen(strcat(lst_files_whole(m).name,'.hdr'), 'r'); %Open file for reading -% for k=1:15 -% tline=fgetl(fid); -% if k==5 -% strscantime=textscan(tline,'%s %s %s %s %s'); -% H_M_S=regexp(cell2mat(strscantime{5}),':','split'); -% frame_start_time(m,1)=str2double(H_M_S(1))*3600+str2double(H_M_S(2))*60+str2double(H_M_S(3)); -% else -% end -% clear tline -% end -% fclose(fid); -% clear k -% end -% clear m -for m=1:size(lst_files_whole,1) - if m==1 - fid = fopen(strcat(lst_files_whole(m).name,'.hdr'), 'r'); %Open file for reading - for k=1:15 - tline=fgetl(fid); - if k==5 - strscantime=textscan(tline,'%s %s %s %s %s'); - H_M_S=regexp(cell2mat(strscantime{5}),':','split'); - world_clock=str2double(H_M_S(1))*3600+str2double(H_M_S(2))*60+str2double(H_M_S(3)); - else - end - clear tline - end - fclose(fid); - end - [duration,ids]=calc_lst_duration(lst_files_whole(m).name,calc_t_rez); - frame_start_time(m,1)=ids(1); - frame_stop_time(m,1)=ids(2); -end -tdiff=frame_start_time(1,1)-world_clock; -frame_start_time=floor(frame_start_time-tdiff); -frame_stop_time=floor(frame_stop_time-tdiff); -acquisition_start=frame_start_time(1,1); - -%%Set start time and stop time relative to Series Start Time -frame_start_time_true=frame_start_time; -frame_stop_time_true=frame_stop_time; -frame_start_time=frame_start_time_true-frame_start_time_true(1,1); -frame_stop_time=frame_stop_time_true-frame_start_time_true(1,1); - -%%Build Time Line (1 if data was collect during that second, 0 if no data -%%was present) -Scanner_Duration((frame_start_time(1)+1):(frame_stop_time(size(frame_stop_time,1))))=0; -for m=1:(size(frame_start_time,1)) - Scanner_Duration((frame_start_time(m)+1):(frame_stop_time(m)))=1; -end -clear m -%Remove short breaks due to framing data -Series_Duration=Scanner_Duration; -System_off=find(Series_Duration==0)'; -for m=1:size(System_off,1) - if Scanner_Duration(System_off(m))==0 && Series_Duration(System_off(m))==0 - if Series_Duration(System_off(m)-1)==1 && Series_Duration(System_off(m)+10)==1 - Series_Duration(System_off(m):System_off(m)+10)=1; - end - else - end -end -% for m=1:size(System_off,1) -% if sum(Series_Duration((System_off(m)-5):(System_off(m))))~=0 && sum(Series_Duration((System_off(m)):(System_off(m)+5)))~=0 -% Series_Duration(System_off(m))=1; -% else -% end -% end -clear m - -end -%important variables -%%Scanner_Duration Series_Duration frame_start_time frame_stop_time -%%lst_files lst_files_whole \ No newline at end of file diff --git a/gamma_bomb_2/code/calc_lst_duration.m b/gamma_bomb_2/code/calc_lst_duration.m deleted file mode 100644 index 3ffa625..0000000 --- a/gamma_bomb_2/code/calc_lst_duration.m +++ /dev/null @@ -1,101 +0,0 @@ -%six bytes ber given event -function [total_time,ids,last_gc]=calc_lst_duration(lst_ffile,t_mks) - -fid = fopen(lst_ffile, 'r', 'l'); %Open file for reading - -if nargin<2, t_mks=1000;end - -%% Find first x-time marks -index=0; -test={}; -t=0; -while 1 - fseek(fid, 6*index, 'bof'); - %start_event=reshape(flipud(dec2hex(fread(fid,6,'uint8')))',1,12); - start_event=reshape(flipud(dec2hex(fread(fid,6,'uint8'),2))',1,12); - if (strcmp(start_event(1,2),'A') && strcmp(start_event(1,3),'0')) || feof(fid) - t=t+1; - test{end+1,1}=start_event; - if t==t_mks || feof(fid) - break - else - index=index+1; - end - else - index=index+1; - end -end - -time={}; -for m=1:size(test,1) - time{end+1,1}=hex2dec(test{m,1}(1,4:end)); -end -errtest=[1;cell2mat(time(2:end))-cell2mat(time(1:end-1))]; -if ~isempty(find(errtest>1)) - start_event=test{test{1}}; - stop_event=test{find(errtest>1)-1}; - start_time=hex2dec(start_event(1,4:end)); - stop_time=hex2dec(stop_event(1,4:end)); - extra_time=floor((stop_time-start_time)*200*10^(-6)); - - start_event=test{find(errtest>1)}; -else - start_event=test{1}; - extra_time=0; -end - - -%% Find last x-time marks -index=1; -test2={}; -t=0; -while 1 - fseek(fid, -6*index, 'eof'); - %stop_event=reshape(flipud(dec2hex(fread(fid,6,'uint8')))',1,12); - stop_event=reshape(flipud(dec2hex(fread(fid,6,'uint8'),2))',1,12); - if index==1 %get last_gc - last_gc=stop_event(1,1); - end - if (strcmp(stop_event(1,2),'A') && strcmp(stop_event(1,3),'0')) || ~ftell(fid) - t=t+1; - test2{end+1,1}=stop_event; - if t==t_mks || ~ftell(fid) - break - else - index=index+1; - end - else - index=index+1; - end -end - -% %% Error correction for frame end (added 05/06/2013) -% time={}; -% for m=1:size(test2,1) -% time{end+1,1}=hex2dec(test2{m,1}(1,4:end)); -% end -% errtest2=[1;cell2mat(time(2:end))-cell2mat(time(1:end-1))]; -% if ~isempty(find(errtest2~=-1)) -% estart_event=test2{test2{1}}; -% estop_event=test2{find(errtest2~=-1,1,'last')-1}; -% estart_time=hex2dec(estart_event(1,4:end)); -% estop_time=hex2dec(estop_event(1,4:end)); -% extra_time=extra_time+floor((estart_time-estop_time)*200*10^(-6)); -% -% stop_event=test2{find(errtest2~=-1,1,'last')}; -% else -% stop_event=test2{1}; -% extra_time=extra_time; -% end - -stop_event=test2{1}; - -start_time=hex2dec(start_event(1,4:end)); -stop_time=hex2dec(stop_event(1,4:end)); - -total_time=floor((stop_time-start_time)*200*10^(-6))+extra_time; - -fclose(fid); - -ids=[start_time;stop_time]*200*10^(-6); -end \ No newline at end of file diff --git a/gamma_bomb_2/code/check_decreasing.m b/gamma_bomb_2/code/check_decreasing.m deleted file mode 100644 index d76e516..0000000 --- a/gamma_bomb_2/code/check_decreasing.m +++ /dev/null @@ -1,84 +0,0 @@ -function [num_deleted, adjusted_data] = check_decreasing(data) -% Check that data is monatonically decreasing after initial hill -% If data prematurely approaches zero, omit these points from analysis - -tail_avg = 0; -temp = zeros(0,2); % Temporary matrix -cutoff_min = 60; % Program starts looking for bad points here -cutoff_max = 200; - -% Calculate the average of all measurements taken after 60 seconds -for m=1:size(data) - if data(m,1) > cutoff_min - temp(end+1,:) = data(m,:); - tail_avg = mean(temp); - end -end -clear temp; -approx_zero = min(tail_avg(1,2)./2, data(end,2)); % Cutoff for what is considered to be "near zero" - -% Add a third column to the data table (1 = keep the point, 0 = delete it) -% Default: keep all points -for m=1:size(data) - data(m,3) = 1; -end - -% Determine the minima within the range of analysis (between time at cutoff_min & -% time at cutoff_max) -for m=2:size(data)-1 - if(data(m,2)cutoff_max) - data(m,3) = 1; - end - if(data(m,2)>approx_zero) - data(m,3) = 1; - end -end - -% Delete points adjacent to the minima found above when necessary -larger = 0; -for m=2:size(data)-1 - if(data(m,3) == 0) - larger = max(data(m+1,2), data(m-1,2)); - end -end -for m=2:size(data)-1 - for m=2:size(data)-1 - if data(m,2) < larger && (data(m+1,3)==0 || data(m-1,3)==0) - data(m,3) = 0; - end - end -end - - -% Store points we'll keep for further analyses in 'adjusted_data' -% Store all deleted points in deleted_points -adjusted_data = zeros(0,2); -deleted_points = zeros(0,2); -for m=1:size(data) - if data(m,3) == 1 - adjusted_data(end+1,:) = data(m,1:2); - else - deleted_points(end+1,:) = data(m,1:2); - end -end -%adjusted_data(end+1,:) = data(end,1:2); % Add final data point to adjusted_data -[num_deleted,~] = size(deleted_points); - -%% Plot to Debug -% figure -% hold on; -% plot(data(:,1),data(:,2), '.r'); -% set(gca, 'XScale', 'log'); -% title('Before'); -% plot(adjusted_data(:,1),adjusted_data(:,2), '.'); -% set(gca, 'XScale', 'log'); -% refline(0, approx_zero); -% refline(Inf, cutoff_min); -% refline(Inf, cutoff_max); -% title('Adjusting for Points that Prematurely Approach Zero'); -% hold off; - -end \ No newline at end of file diff --git a/gamma_bomb_2/code/convertshc.m b/gamma_bomb_2/code/convertshc.m deleted file mode 100644 index 1e5fa8a..0000000 --- a/gamma_bomb_2/code/convertshc.m +++ /dev/null @@ -1,40 +0,0 @@ -function convertshc(hc_file,hp_hc) - -fid=fopen(hc_file); -t_hc_line=[]; -while 1 - hc_line=fgetl(fid); - if hc_line==-1,break;end %EOF - if strcmpi(hc_line(1:10), 'Head Curve') - t_hc_line(end+1,:)=sscanf(hc_line,'Head Curve %i %i %i %i %i %i %i %i %i %i %i ET=%f gcerrors=%i')'; - end -end -fclose(fid); - - -if nargin>1 - t_hc_line2=csvread(hp_hc, 1, 0); - %%make sure HistPET is always increasing - for m=2:size(t_hc_line2,1); - if t_hc_line2(m,1)size(t_hc_line2,1) %% HistPET version has fewer times...lets just assume they increase by 1 sec as one might assume. - offset=size(t_hc_line,1)-size(t_hc_line2,1); - t_hc_line2=[t_hc_line2;[[1:offset]'*1000,zeros(offset,12)]]; - end - t_hc_line(:,1)=t_hc_line2(1:size(t_hc_line,1),1); -end - - -if ~isempty(t_hc_line) - fid1=fopen(hc_file, 'wt'); - fprintf(fid1,'Time(msec),CFD,StoredXY,QualEvent,RIO_Singles,Delays,RIO_Delay,Prompts,RIO_Prompt,ERS_Prompt,ERS_Delays,ERS_Singles,SyncErrCnt\n\r'); - for m=1:size(t_hc_line,1) - fprintf(fid1,[sprintf('%i,%i,%i,%i,%i,%i,%i,%i,%i,%i,%i,%i,%i',t_hc_line(m,:)) '\n\r']); - end - fclose(fid1); -end - -end \ No newline at end of file diff --git a/gamma_bomb_2/code/dos_ex.m b/gamma_bomb_2/code/dos_ex.m deleted file mode 100644 index b315f1f..0000000 --- a/gamma_bomb_2/code/dos_ex.m +++ /dev/null @@ -1,9 +0,0 @@ -function dos_ex(dos_string) -%%I'm not sure why, but the dos/system function crashes when it is using -%%lmsort, but seems to be okay with !. Let's make an easy function to do -%%this - - -evalin('base',[strcat('!',dos_string)]) - -end \ No newline at end of file diff --git a/gamma_bomb_2/code/feng_fit3.m b/gamma_bomb_2/code/feng_fit3.m deleted file mode 100644 index bee139d..0000000 --- a/gamma_bomb_2/code/feng_fit3.m +++ /dev/null @@ -1,246 +0,0 @@ -function [fitobject, fit_metric, models]=feng_fit3(bld_ffile,endtime,met_fit,exclude_from,sf,options) -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%% SCRIPT TO FIT BLOOD DATA WITH FENG & L3 MODELS -%% Written by Daniel Chonde -%% Req:blood file (*.bld) -%% Output: blood file (*_feng.bld) -%% Written: 08/12/2012 Updated: -% Will also now return statistical gof data (goodness of fit) -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% This script reads in a .bld file and fits Feng's model to the blood data. -% -% example: -% feng_fit(bld_ffile,endtime) -% bld_ffile=path to a .bld file -% endtime=time to interpolate model to in seconds -% sf=number of decimal places to write to file -% -% feng_fit(bld_ffile) fits the model and interpolates the data to the -% last timepoint in the bld file -% -% feng_fit(bld_ffile,endtime) fit the model and interpolate the data to -% an arbitrary timepoint given by endtime, where endtime is in seconds. -% -fit_metric = zeros(0,1); - -if nargin<4, exclude_from='';end - -A=dlmread(bld_ffile,'\t',1,0); -if A(1,1)~=0, A=[0 0;A];end %add t=0 activity=0, if not already present - -if exclude_from - A(find(A(:,1)>exclude_from,1,'first'):end,:)=[]; -end - -if nargin<5 || isempty(sf) %%Figure out how many decmal places are used - %fid=fopen(bld_ffile,'r') - %test=textscan(fid,'%s') - %fclose(fid) - teststr=num2str(A(end,2)); - splits=regexp(teststr,'\.','split'); - sf=size(splits{2},2); -end -names=textread(bld_ffile,'%s',1,'delimiter','\n'); -hdrnames=regexp(names,'\t','split'); - -time=A(:,1);plasma=A(:,2); -num_points = length(time); -assignin('base','num_points',num_points); -if nargin<2, endtime=max(time);end -if nargin<3, met_fit=fit([1:10]',ones(1,10)','poly1');end % if not met data, construct a straight line -%% Remove points after 5 minutes that go up greater than 5% -cutofft=find(time>300,1,'first'); -iter=1; -num_outliers = 0; -while 1 - kvals=[zeros(cutofft-1,1);(plasma(cutofft:end).*met_fit(time(cutofft:end))-plasma(cutofft-1:end-1).*met_fit(time(cutofft-1:end-1)))./(plasma(cutofft-1:end-1).*met_fit(time(cutofft-1:end-1)))>0.05]; - if iter==1 - ktime=time(find(kvals)); - kplasma=plasma(find(kvals)); - else - ktime=[ktime;time(find(kvals))]; - kplasma=[kplasma;plasma(find(kvals))]; - end - num_outliers = length(kvals); - time(find(kvals))=[]; - plasma(find(kvals))=[]; - iter=iter+1; - if ~any(kvals) - num_outliers = 0; - break; - end -end -assignin('base', 'num_outliers', num_outliers); - -%%Quick Linear interpolation so we can have that file (used by Jean Logan) -plasma_linear=interp1(time,plasma,time(1):endtime,'linear','extrap'); -plasma_linear(plasma_linear<0)=0; - - -%% Fit Models -models={'feng_model_aif(a,b,c,d,e,f,g,x)','Feng' - 'lin_3exp_model_aif(a,b,c,d,e,f,g,h,jj,x)','L3Exp'}; - -%Check if there is an options.Choice--this is my workaround to allow you to -%select a single fit model -if exist('options','var') && isfield(options,'Choice'),models(:,~cellfun(@(x) strcmpi(options.Choice,x),models(:,2)))=[];end - -for ppp=1:size(models,1) - ffun=fittype(models{ppp,1}); - if nargin<6 - options=fitoptions(ffun); - options.MaxIter=1000; - options.MaxFunEvals=1000; - switch models{ppp,2} - case 'Feng' - options.Lower=[0 -1 0 -1 0 -1 0]; - options.Upper=[500 0 500 0 500 0 time(find(plasma==max(plasma)))+30]; - %options.Weights=sqrt(plasma)/sum(sqrt(plasma)); - case 'L3Exp' - options.MaxFunEvals=1000;options.Lower=[0 0 0 0 0 0 0 0 0 ]; - options.Upper=[1000 1 1000 1 1000 1 30 time(find(plasma==max(plasma))) 30]; - end - end - - for m=1:2 - if m==1 - switch models{ppp,2} - case 'Feng' - if isempty(options.StartPoint),options.StartPoint=[251/60,-2.1/60,5.8,-0.0104/60,1.879,-0.219/60,time(find(plasma==max(plasma)))-20];end - case 'L3Exp' - if isempty(options.StartPoint),object.StartPoint=[60 .1 90 11 138 .02 10 1 1];end - end - else - switch models{ppp,2} - case 'Feng' - options.StartPoint=[fitobject.a,fitobject.b,fitobject.c,fitobject.d,fitobject.e,fitobject.f,fitobject.g]; - case 'L3Exp' - options.StartPoint=[fitobject.a,fitobject.b,fitobject.c,fitobject.d,fitobject.e,fitobject.f,fitobject.g,fitobject.h,fitobject.jj]; - end - end - - [fitobject,gof]=fit(time,plasma,ffun,options); - display(['Iteration ' num2str(m)]) - display(gof) - display(fitobject) - end - - fplasma=fitobject(time(1):endtime); - - - - % ***SUBPLOT #1*** Plot data on Log Scale x-axis - assignin('base', 'kplasma', kplasma); - figure; subplot(2,2,1); - if(size(kplasma)>0) - plot(time,plasma,'xb');hold on;plot(ktime,kplasma,'*g');fplot(fitobject,[time(1) endtime],0.0001,'r');set(gca,'XScale','log'); - legend('Data','Outliers','Fitted Curve'); - legend off; - else - plot(time,plasma,'xb');hold on;plot(ktime,kplasma,'*g');plot(fitobject);title(strcat(models{ppp,2},' Model Fit')); - end - title(strcat(models{ppp,2},' Model Fit (Log Scale)')); - xlabel('Time (s)','FontSize',7);ylabel('kBq/ml','FontSize',7); - - - - % ***SUBPLOT #2*** Plot data on regular x-axis - subplot(2,2,2); - if ~isempty(kplasma) - plot(time,plasma,'xb');hold on;plot(ktime,kplasma,'*g');plot(fitobject);title(strcat(models{ppp,2},' Model Fit')) - legend('Data','Outliers','Fitted Curve'); - else - plot(time,plasma,'xb');hold on;plot(fitobject); - legend('Data', 'Fitted Curve', 'Location', 'best'); - end - xlabel('Time (s)','FontSize',7);ylabel('kBq/ml','FontSize',7); - - - - % Calculate integral for raw data (Crude approximation on graph) - int_of_raw = zeros(length(time),1); - raw_data = cat(2, time, plasma); - for m=1:length(int_of_raw) - if m==1 - int_of_raw(m,1) = (raw_data(m,1)).*(raw_data(m,2)); - else - int_of_raw(m,1) = ((raw_data(m,1)-raw_data(m-1,1)).*(raw_data(m,2))) + int_of_raw(m-1,1); - end - assignin('base', 'int_of_raw', int_of_raw); -% int_of_raw(m,1) = trapz(raw_data(1:m,2)); - end - % Calculate integral for fitted data - fit_data = zeros(length(fplasma),2); - for m=0:length(fplasma)-1 - fit_data(m+1,1) = m; - fit_data(m+1,2) = fplasma(m+1,1); - end - int_of_fit = zeros(length(fplasma),1); - for m=1:length(int_of_fit) - int_of_fit(m,1) = trapz(fit_data(1:m,2)); - end - - - % ***SUBPLOT #3*** Plot the integral of the fit compared to the - % integral of the raw data - subplot(2,2,3); - hold on; - plot(fit_data(:,1),int_of_fit(:,1),'r'); - plot(time,int_of_raw,'b'); - % Adjust Plot Format - title('Integral of Fit vs. Crude Approximation', 'FontSize', 9, 'FontWeight', 'bold'); - %legend('Integral of Fit', 'Crude Approximation', 'Location', 'best'); - xlabel('Time (s)', 'FontSize', 7); ylabel('Integral', 'FontSize', 7); - set(gca, 'FontSize', 7, 'XScale', 'log'); - - - - % ***SUBPLOT #4*** Integrals on linear x-scale - subplot(2,2,4); - plot(fit_data(:,1),int_of_fit(:,1),'r'); - hold on; - plot(time,int_of_raw,'b'); - title('Integral of Fit vs. Crude Approximation', 'FontSize', 9, 'FontWeight', 'bold'); - legend('Integral of Fit', 'Crude Approximation', 'Location', 'SouthEast'); - xlabel('Time (s)', 'FontSize', 7); ylabel('Integral', 'FontSize', 7); - set(gca, 'FontSize', 7); - - fit_metric(end+1) = evaluate_fit(int_of_fit,int_of_raw); - - - - - % Create textbox to display "fit metric" - annotation(gcf,'textbox',... - [0.841619718309864 0.00418410041840576 0.147112676056338 0.0719487425300029],... - 'String',{'Fit Metric: ' num2str(fit_metric(1, ppp))},... - 'FontWeight','bold',... - 'FitBoxToText','on',... - 'BackgroundColor',[1 1 1]); - - % Save the graph as an image - drawnow; - set(gcf, 'Position', [50 50 840 630]); - saveas(gcf,fullfile(pwd,regexprep(bld_ffile,'\.bld',['_' models{ppp,2} '.bmp']))) - - % Save the BLD files - [bldpath,bldname,bldext]=fileparts(bld_ffile); - print_list(hdrnames{:},[[time(1):endtime]', fplasma],fullfile(bldpath,strcat(bldname,['_' models{ppp,2}],bldext)),sf); -end - print_list(hdrnames{:},[[time(1):endtime]', plasma_linear'],fullfile(bldpath,strcat(bldname,'_linear',bldext)),sf); - disp(fit_metric); -end - - - -function fit_metric = evaluate_fit(int_of_fit,int_of_raw) -fit_metric = nan; -% Compares the integral of the fit data and the integral of the raw data to -% determine goodness of fit. -% Returns a "fit metric" - -% Calculate integral for the fit data -fit_metric = log10(abs(int_of_raw(end,1) - int_of_fit(end,1))); -disp 'FIT METRIC'; -disp(fit_metric); -end \ No newline at end of file diff --git a/gamma_bomb_2/code/feng_model_aif.m b/gamma_bomb_2/code/feng_model_aif.m deleted file mode 100644 index 7651984..0000000 --- a/gamma_bomb_2/code/feng_model_aif.m +++ /dev/null @@ -1,22 +0,0 @@ -function aif=feng_model_aif(A1,L1,A2,L2,A3,L3,tau,t) -%%Model of AIF bolus based on Feng et al. Int J Biomed Comput (1993) -% Y(t)=0, t=tau -% -% [A1]=uCi/ml/sec -% [A2]=[A3]=uCi/ml -% [t]=sec -% [tau]=sec -% [L1]=[L1]=[L1]=1/sec - -aif=zeros(size(t)); - -pre= t=tau; - -aif(pre)=0; -aif(post)=(A1*(t(post)-tau)-A2-A3).*exp(L1*(t(post)-tau))+... - A2*exp(L2*(t(post)-tau))+A3*exp(L3*(t(post)-tau)); - -end \ No newline at end of file diff --git a/gamma_bomb_2/code/final_report.m b/gamma_bomb_2/code/final_report.m deleted file mode 100644 index 07d61c5..0000000 --- a/gamma_bomb_2/code/final_report.m +++ /dev/null @@ -1,180 +0,0 @@ -%% Gamma Bomb! 2.0 Blood Data Report -%% Participant Information -% This section displays important demographic information about the -% participant and gives some generic information about the present study. -%% -[~,~,chart] = xlsread(dose_info); -assignin('base', 'chart', chart); -name = ''; -start=regexp(participant, ligand, 'start'); -for m=start(1,1):(start(1,1)+11) - name = strcat(name, participant(m)); -end -disp(strcat('Participant #:.......',name)); -disp(strcat('Dose Info Sheet:.....',dose_info)); -hour = floor(toi_time./3600); -minute = floor(mod(toi_time,3600) ./ 60); -second = floor(mod(mod(toi_time,3600),60)); -t_hour = num2str(hour); -t_minute = num2str(minute); -t_second = num2str(second); -if hour < 10 - t_hour = strcat('0',num2str(hour)); -end -if minute < 10 - t_minute = strcat('0',num2str(minute)); -end -if second < 10 - t_second = strcat('0',num2str(second)); -end -disp(' '); -disp(strcat('TOI:.........................',t_hour,':',t_minute,':',t_second)); -disp(strcat('BAT Offset...................', num2str(batoffset), ' sec')); -if isnan(chart{16,3}) - radiotracer = tracer; -else - radiotracer = chart{16,3}; -end -disp(['Radioligand:.................' radiotracer '- ' ligand]); -disp(' '); -if isnan(chart{2,10}) - weight = 'Unavailable'; -else - weight = [num2str(chart{2,10}) ' lbs']; -end -disp(strcat('Participant Weight...........', weight)); -if isnan(chart{5,10}) - height = 'Unavailable'; -else - height = [num2str(chart{5,10}) ' in']; -end -disp(strcat('Participant Height...........', height)); -if isnan(chart{8,10}) - sex = 'Unavailable'; -else - sex = upper(chart{8,10}); -end -disp(strcat('Participant Sex..............', sex)); -if isnan(chart{11,10}) - byear = 'Unavailable'; -else - byear = num2str(chart{11,10}); -end -disp(strcat('Participant Birth Year.......', byear)); -disp(' '); - -best = min(abs(models_gof)); -best_model_i = nan; -for m=1:length(models_gof) - if models_gof(m) == best - best_model_i = m; - end -end -best_model = models(best_model_i, 2); -disp(strcat('Best TAC Model...............', best_model)); -best_fit_i = nan; -best = max(abs(GOF_TABLE(:,1))); -for m=1:length(GOF_TABLE) - if GOF_TABLE(m) == best - best_fit_i = m; - end -end -best_fit = upper(pf_fits(best_fit_i,1)); -disp(strcat('Best Parent Fraction Fit.....', best_fit)); -%% -% -------------------------------------------------------------------------------------------------------------------------- -%% Blood Time Activity Curve (TAC) Fits -% _Current analysis uses two kinetic models to interpolate the Time Activity -% Curve:_ -% -% # Feng Fit -% # Linear to 3-Exponential (L3Exp) Fit -% -% _During this analysis, outliers were removed prior to curve-fitting. Outliers are -% classified as any point after 5 minutes that increased greater than 5 -% percent OR as any point that prematurely approaches zero._ -% -%% -disp(['Number of Outliers: ' num2str(num_outliers+num_deleted)]); -disp(['Number of Points Used in Analysis: ' num2str(num_points)]); -for m=1:3 - disp(' '); -end -%% -% -% *Goodness of Fit* -% -% The table below summarizes the goodness of fit of the two TAC Models. -% The fit metric represents the space between the fitted curve and the raw -% data. It is calculated according to the equation: -% -%% -% -% $M = \log ( |\int_{0}^{t} x dx - \int_{0}^{t} x' dx |)$ -% -% _where:_ -% -% * $M$ = _goodness of fit_ -% * $t$ = _duration of the blood data (s)_ -% * $x$ = _raw data points_ -% * $x'$ = _fitted data_ -% -% _According to the metric, a perfect fit would yeild M = 0. The "best fit" -% will have a value M closest to 0._ -% -%% -hdrs = ''; -types = models(:,2); -for m=1:size(types) - if m==1 - hdrs = [types(1)]; - else - hdrs = [hdrs ' ' types(m)]; - end -end -printmat(models_gof, 'Model Goodness of Fit', 'Fit_Metric', 'Feng_Fit L3Exp_Fit'); -disp(' '); -disp(['Best TAC Model: ' best_model]); - -%% -% *Feng Model* -% -% <> -% -% *L3Exp Model* -% -% <> -% -%% -% -------------------------------------------------------------------------------------------------------------------------- -%% Summary of Parent Fraction Fit Statistics -% _The Parent Fraction has been fitted to four different curves:_ -% Hill Fit, Linear to Exponential Fit (Lin2Exp), Exponential Fit (Exp), and -% Power Fit. -% -%% -printmat(GOF_TABLE, 'Parent Fraction Fit Stats', 'Hill_Fit Lin2Exp_Fit Exponential_Fit Power_Fit',... - 'Rsquare Adj_RSquare DFE'); -disp(' '); -disp(['Best Parent Fraction Fit: ' best_fit]); -cd (csv_path); -%% Parent Fraction Fits -% -% <> -% -%% Histograms for Parent Fraction Data Points -% _Each histogram below represents the blood data for single the Parent -% Fraction Data Point indicated by the graph's title. Red bars show -% metabolite data and blue bars show the parent data._ -% -% <> -% -% <> -% -%% -if population_bool - disp('Parent Fraction calculated from Population data. Histograms unavailable.'); -end - -%% -% GammaBomb! 2.0 - Developed at the Hooker Research Group - MGH Martinos Center - Boston, MA \ No newline at end of file diff --git a/gamma_bomb_2/code/find_BAT.m b/gamma_bomb_2/code/find_BAT.m deleted file mode 100644 index 127388e..0000000 --- a/gamma_bomb_2/code/find_BAT.m +++ /dev/null @@ -1,89 +0,0 @@ -function InjectionStart=find_BAT(head_curve1,crop_counts,save_img) -if nargin<2 - %crop_counts=300; - crop_counts=300; -end -if nargin<3,save_img=1;end -if length(head_curve1)>crop_counts -head_curve=head_curve1(1:crop_counts); -else - head_curve=head_curve1; -end -%% Clean up the head curve so we can estimate the injection time -signoise=nanmean(head_curve(1:5)); -if isnan(signoise), signoise=0;end -head_curve=head_curve-signoise; %consider the first 5 points as noise -head_curve(head_curve<0)=NaN; -head_curve=pchip(1:size(head_curve,1),head_curve,[1:size(head_curve,1)]'); -%Quick moving average filter to remove the noise -b=ones(1,16)/16; -a=1; -filteredData=filter(b,a,head_curve); -dtest=(filteredData(2:end)-filteredData(1:end-1))./(filteredData(2:end)+filteredData(1:end-1)); -fdtest=filter(ones(1,5)/5,a,dtest); -sfdtest=fdtest;sfdtest(filteredData<1000)=0; -[~,I]=max(sfdtest); -assignin('base', 'a', a); -assignin('base', 'b', b); -%g=@(a,b,c,d,e,x) model_aif(a,b,c,d,e,x); -%ffun=fittype(g); -ffun=fittype('model_aif(a,b,c,d,e,x)'); - -assignin('base', 'a_after', a); -assignin('base', 'ffun', ffun); - -options=fitoptions(ffun); - -assignin('base', 'options', options); - -%options.StartPoint=[1,30,min(head_curve(head_curve>0)),(max(head_curve(head_curve>0 &head_curve0)))/10,(head_curve(ceil(length(head_curve)*1/2))-head_curve(end))/length(head_curve)]; -x_val=[0:size(head_curve,1)-1]'; -x_val_fit=x_val(head_curve>=0 &head_curve=0 &head_curve1000 -%options.Lower=[x_val_fit(find(head_curve_fit>=100,1,'first'))-.1 0 min(head_curve_fit) -(max(head_curve_fit)-min(head_curve_fit)) -(max(head_curve_fit)-min(head_curve_fit))]; -%options.Upper=[x_val_fit(find(head_curve_fit>=100,1,'first')) 300 max(head_curve_fit) (max(head_curve_fit)-min(head_curve_fit)) (max(head_curve_fit)-min(head_curve_fit))]; - -%% -[fitobject,gof]=fit(x_val_fit,head_curve_fit,ffun,options); - -assignin('base', 'head_curve_fit', head_curve_fit); -assignin('base', 'x_val_fit', x_val_fit); -assignin('base', 'fitobject', fitobject); -assignin('base', 'a_fitobject', fitobject.a); - -if fitobject.a>0 && fitobject.d>0 -InjectionStart=ceil(fitobject.a); -else -InjectionStart=0; -end -if save_img -figure('Position',[680 558 560*2 420]);%figure; -for gg=1:2 -subplot(1,2,gg);plot(0:length(head_curve)-1,head_curve); -hold on -fplot(fitobject,[0,length(head_curve)-1],'r'); -set(get(gca,'XLabel'),'String','Time (sec)') -set(get(gca,'YLabel'),'String','Counts') -set(get(gca,'Title'),'String',['Bolus Arrival Time Pick-off T=' num2str(InjectionStart)]) -if fitobject.a>0 && fitobject.d>0 && gg==1 -xlim([0,300]) -else - xlim([0,20+fitobject.a]) -end -end -legend('head curve','fit','Location','NorthWest'); -%saveas(gcf,'BAT.pdf'); -set(gcf,'PaperPositionMode','auto','PaperOrientation','landscape');print(gcf,'-dpdf','BAT.pdf') -end - -end \ No newline at end of file diff --git a/gamma_bomb_2/code/gamma_bomb.m b/gamma_bomb_2/code/gamma_bomb.m deleted file mode 100644 index 35e1c04..0000000 --- a/gamma_bomb_2/code/gamma_bomb.m +++ /dev/null @@ -1,605 +0,0 @@ -% patient_AIF_gamma_read.m -% Author: Spencer Bowen -% Lab Group: Catana -% Date: 09-06-12 -% Read automatic gamma counter data for AIF arterial and venous whole blood and plasma samples and -% sort into a readable format. -% inputs: template_file = name of csv file with automatic gamma counter sample information -% csv_file = name, including directory, of the csv file with the first set of -% sample measurements (assumes csv files numeric with sequential values) -% csv_out = name of the csv file to write the compiled data to -% -% -% Modified By Tom Morin -% 5/01/15 - Modified to print results to a single html/PDF file -% 5/20/15 - Added Power, Exponential, & Hill Fits to Parent Fraciton Data -% 5/21/15 - Now displays BAT Offset (Bolus Arrival Time) in final report -% and accounts for BAT in .bld output files -% 5/22/15 - More aesthetically pleasing final reports -% 5/27/15 - Delete points that prematurely approach zero in TAC -% 6/2/15 - Added histograms to show origin of parent fraction data -% 6/3/15 - Display integrals with Fit/L3Exp Fits & calculate custom metric -% 6/4/15 - State which fits/models are best in the final report. - -function gamma_bomb(template_file, csv_path, tracer, met_correct, dose_info, report_format, population_pf, bk_cps_value) -%% Check to see if a Final Report has already been created -exit_gamma_bomb = 0; -exit_gamma_bomb = check_overwrite(csv_path, population_pf); -if(exit_gamma_bomb == 1) - disp('Exiting Gamma Bomb 2.0'); - return; -end - -assignin('base', 'population_bool', population_pf); -assignin('base', 'tracer', tracer); - -%% Define Constants -switch lower(tracer) - case {'f-18'} %F-18 - half_life=60*120; - case {'c-11'} %c-11 - half_life=20.38*60; - otherwise - error('unknown tracer. Contact Ciprian') -end - -%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% define constants -% number of positions in the gamma counter racks -pos_num = 10; -max_samp= 10000; -csv_dat_str = '%d%s%s%d%d%d%d%d%f%s%f%f%f%s'; -pa_col = 1; -wa_col = 2; -pv_col = 3; -wv_col = 4; -out_str = '%d,%d,%d,%s,%f,%d'; -num_out_col = 6; -% strings describing output columns -out_hdr_str = {'Run ID','Rack','Pos','Measurement Time',... - 'Duration (sec)','F-18 350-600 keV Counts'}; - -%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% load parameters -current_dir = pwd;cd(csv_path); - -%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%% Read in Template File -h = waitbar(0, 'Reading in Template Files - 0%'); -[~,tube_code_mx,~]=xlsread(template_file,1); -[col_time_mx,~,~]=xlsread(template_file,2); -[tube_vol_mx,~,~]=xlsread(template_file,3); -[TOI,~,~]=xlsread(template_file,4); -toi_vec=datevec(TOI); -assignin('base', 'TOI', TOI); -toi_time=toi_vec(4)*3600+toi_vec(5)*60+toi_vec(6); -%assignin('base', 'toi_vec', toi_vec); -assignin('base', 'toi_time', toi_time); -InjectionStart=0; -update_waitbar(0.05, h, 'Calculating BAT (This may take a few minutes.)'); -if exist(fullfile(home_dir,'PET','Frame1.lst.hdr'),'file') - if exist(fullfile(home_dir,'PET','Sum','TAC_workspace.mat'),'file') - display('using TAC to determine the start time and align the data to the PET') - load(fullfile(home_dir,'PET','Sum','TAC_workspace.mat'),'InjectionStart'); - load(fullfile(home_dir,'PET','Sum','TAC_workspace.mat'),'Series_Duration'); - else - display('No TAC workspace found--Reconstrucing Frame1.lst to get injection Start') - c_dir=pwd; - cd(fullfile(home_dir,'PET')); - [Scanner_Duration, Series_Duration, frame_start_time, ... - frame_stop_time, acquisition_start, lst_files, ... - lst_files_whole]=TAC_sort_lst_mode(fullfile(home_dir,'PET')); - TAC_process_data('Frame1'); - update_waitbar(0.1, h, 'Generating Head Curve...'); - file=fullfile(pwd, strcat('Frame1','_head_curve.s.hc')); %build fullfile name - [head_curve,head_curve_trues]=generate_head_curve(file,Series_Duration); - update_waitbar(0.2, h, 'Determining Injection Start Time'); - - assignin('base', 'head_curve_trues', head_curve_trues); - - InjectionStart=find_BAT(head_curve_trues); - update_waitbar(0.22, h, 'Deleting Frames'); - delete('Frame1*.lor*','Frame1*.s*') - update_waitbar(0.24, h,'Changing Directory'); - cd(c_dir) - end - update_waitbar(0.26, h, 'Dipslaying Injection Start & Series Duration'); - display(['Injection Start: ' num2str(InjectionStart)]); - display(['Series Duration: ' num2str(length(Series_Duration))]); - update_waitbar(0.28, h, 'Determining extrap_time'); - extrap_time=length(Series_Duration); - update_waitbar(0.3, h, 'Crunching some more numbers...'); - lst_hdr=read_i_hdr(fullfile(home_dir,'PET','Frame1.lst.hdr')); - lst_startvec=datevec(lst_hdr.StudyTimeHhMmSs); - lst_start=lst_startvec(4)*3600+lst_startvec(5)*60+lst_startvec(6); - batoffset=(lst_start+InjectionStart)-toi_time; - toi_time=toi_time+batoffset; - assignin('base','InjectionStart', InjectionStart); - assignin('base', 'lst_start', lst_start); - assignin('base', 'batoffset', batoffset); - assignin('base', 'toi_time', toi_time); -end -col_time_mx(1,:)=[]; -tube_vol_mx(1,:)=[]; -tot_samps=size(unique(find(~strcmpi(tube_code_mx,'E'))),1); %total number of non E -datev=datevec(col_time_mx); -secv=datev(1:end,4)*3600+datev(1:end,5)*60+datev(1:end,6); -col_time_mx=reshape(secv,size(col_time_mx))-toi_time; -%% - -csv_files=ls(fullfile(csv_path,'*.csv')); -[~,csv_name,~]=fileparts(csv_files(1,:)); -%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% open and read in csv files (filenames should end in consecutive numbers) -update_waitbar(0.35, h, 'Reading in CSV files'); -csv_tot = 0; -csv_start_num = str2num(csv_name); -csv_samp_tot = 0; -fail=0; -readon = 1; -while readon - csv_cur= fullfile(csv_path,[num2str(csv_start_num+csv_tot,'%06.0f') '.csv']); - - if ~exist(csv_cur,'file'); - fail=fail+1; - if fail==2; - break - end - else - - fid_in = fopen(csv_cur,'r'); - csv_dat= textscan(fid_in,csv_dat_str,'delimiter',',','HeaderLines',1,'EndOfLine','\n'); - fclose(fid_in); - csv_samp_tot = csv_samp_tot + length(csv_dat{1}); - if csv_tot==0 - csv_comp_dat = csv_dat; - else - rack_offset=csv_comp_dat{6}(end); - for cnti=1:14 - if cnti==6 - csv_comp_dat{cnti} = [csv_comp_dat{cnti}; csv_dat{cnti}+rack_offset]; - else - csv_comp_dat{cnti} = [csv_comp_dat{cnti}; csv_dat{cnti}]; - end - end - end - end - - if csv_samp_tot>=tot_samps - readon=0; - else - csv_tot= csv_tot+1; - end - -end -assignin('base', 'csv_comp_dat', csv_comp_dat); - -%% Construct other necessary matrices -%(Take csv data and format it to tables that are the same dimensions as the one in the template file) -%Rack=row -%Pos=col -update_waitbar(0.4, h, 'Constructing matrices'); -count_time_mx=zeros(size(tube_code_mx)); -count_dur_mx=zeros(size(tube_code_mx)); -t_counts_mx=zeros(size(tube_code_mx)); -dev_cpm_mx=zeros(size(tube_code_mx)); -t_counts_error_mx=zeros(size(tube_code_mx)); - -for m=1:size(csv_comp_dat{1},1) - dv=datevec(csv_comp_dat{3}(m)); - count_time_mx(csv_comp_dat{6}(m),csv_comp_dat{8}(m))=dv(4)*3600+dv(5)*60+dv(6); - count_dur_mx(csv_comp_dat{6}(m),csv_comp_dat{8}(m))=csv_comp_dat{9}(m); - t_counts_mx(csv_comp_dat{6}(m),csv_comp_dat{8}(m))=csv_comp_dat{11}(m); - dev_cpm_mx(csv_comp_dat{6}(m),csv_comp_dat{8}(m))=csv_comp_dat{12}(m); - t_counts_error_mx(csv_comp_dat{6}(m),csv_comp_dat{8}(m))=csv_comp_dat{13}(m); -end - -%% Construct blank scan matrix -Bpos=find(strcmpi(tube_code_mx','B')); -display(nargin) -if nargin==8 - bk_mx=ones(size(count_time_mx))*bk_cps_value; -elseif ~isempty(Bpos) - blanks=zeros(size(tube_code_mx))'; - start_pos=1; - t_counts_mx_p=(t_counts_mx./count_dur_mx)'; - for m=1:length(Bpos) - blanks(start_pos:Bpos(m))=t_counts_mx_p(Bpos(m)); - start_pos=Bpos(m)+1; - end - %finish rest of positions - blanks(blanks==0)=t_counts_mx_p(Bpos(end)); - bk_mx=blanks'; -else - bk_mx=zeros(size(count_time_mx)); -end - -%% Decay Correct Data -update_waitbar(0.45, h, 'Decaying Data'); -elapsed_t_mx=count_time_mx-toi_time; -cps_dkc_mx=(t_counts_mx./count_dur_mx-bk_mx).*exp(log(2).*elapsed_t_mx/half_life); -kbq_ml_dkc_mx=cps_dkc_mx./tube_vol_mx*60/40; -%% Breakdown for Metabolite Correction -% -final_labels={}; -file_name={}; -if ~met_correct - [tube_struct,t_labels]=sort_tube_codes(tube_code_mx,col_time_mx,kbq_ml_dkc_mx); -else - kbq_ml_dkc_met_mx=kbq_ml_dkc_mx; - kbq_ml_dkc_met_mx(find(strcmpi(tube_code_mx,'M')|strcmpi(tube_code_mx,'P')))=kbq_ml_dkc_mx(strcmpi(tube_code_mx,'M')|strcmpi(tube_code_mx,'P')).*tube_vol_mx(strcmpi(tube_code_mx,'M')|strcmpi(tube_code_mx,'P')); - [tube_struct,t_labels]=sort_tube_codes(tube_code_mx,col_time_mx,kbq_ml_dkc_met_mx); - % Metabolite Correction - u_times=unique(tube_struct.M(:,1)); - for m=1:size(u_times,1) - mets(m,1)=sum(tube_struct.M(find(tube_struct.M(:,1)==u_times(m)),2)); - parent(m,1)=sum(tube_struct.P(find(tube_struct.P(:,1)==u_times(m)),2)); - if isfield(tube_struct,'PA') - plasma(m,1)=tube_struct.PA(find(tube_struct.PA(:,1)==u_times(m)),2); - elseif isfield(tube_struct,'PV') - plasma(m,1)=tube_struct.PV(find(tube_struct.PV(:,1)==u_times(m)),2); - else - error('For Metabolite correction either arterial or venous samples must be drawn. Contact Ciprian'); - end - end - total=mets+parent; - met_ratio=parent./total; - tube_struct.M=[0 1;u_times met_ratio]; - tube_struct=rmfield(tube_struct,'P'); - t_labels(cellfun(@(x) strcmpi(x,'P'),t_labels))=[]; -end - -update_waitbar(0.5, h, 'Generating TAC Data...'); -for m=1:size(t_labels,1) - switch lower(t_labels{m}) - case 'pa' - final_labels=[final_labels,{'sample-time[seconds]';'plasma[kBq/cc]'}]; - file_name=[file_name,{'plasma_art.bld'}]; - case 'pv' - final_labels=[final_labels,{'sample-time[seconds]';'plasma[kBq/cc]'}]; - file_name=[file_name,{'plasma_ven.bld'}]; - case 'wa' - final_labels=[final_labels,{'sample-time[seconds]';'whole-blood[kBq/cc]'}]; - file_name=[file_name,{'whole-blood_art.bld'}]; - case 'wv' - final_labels=[final_labels,{'sample-time[seconds]';'whole-blood[kBq/cc]'}]; - file_name=[file_name,{'whole-blood_ven.bld'}]; - case 'm' - final_labels=[final_labels,{'sample-time[seconds]';'parent-fraction[1/1]'}]; - file_name=[file_name,{'parentfraction.bld'}]; - otherwise - end -end - -if (exist(fullfile(pwd,'parentfraction.bld'),'file') || ~met_correct) && (exist(fullfile(pwd,'whole-blood_ven.bld'),'file') || exist(fullfile(pwd,'whole-blood_art.bld'),'file') || exist(fullfile(pwd,'plasma_ven.bld'),'file') || exist(fullfile(pwd,'plasma_art.bld'),'file')) - qans=questdlg('Parent-Fraction and Blood .bld files present. Overwrite?','Blood GUI','Yes','No','No'); -else - qans='yes'; -end - -if strcmpi(qans,'yes') - % Delete points that approach zero prematurely - [num_deleted, temp] = check_decreasing(tube_struct.PA); - tube_struct.PA = []; - tube_struct.PA = temp; - % Print Preliminary Data to BLD files - for m=1:size(t_labels,1) - if ~(population_pf && strcmp(file_name{m},'parentfraction.bld')) % Check to see if Population Parent Fraction has already been created - print_list(final_labels(:,m), tube_struct.(t_labels{m}), fullfile(pwd,file_name{m}), 4); - figure; - if ~strcmpi(t_labels{m},'m') && tube_struct.(t_labels{m})(1,1)~=0 && tube_struct.(t_labels{m})(1,2)~=0 - tube_struct.(t_labels{m})=[0 0;tube_struct.(t_labels{m})]; - end - - % PCHIP Fit of Parent Fraction Data - p_cfit=fit(tube_struct.(t_labels{m})(:,1),tube_struct.(t_labels{m})(:,2),'pchipinterp'); - subplot(2,1,1); - plot(tube_struct.(t_labels{m})(:,1),tube_struct.(t_labels{m})(:,2),'r*'); - hold on; - fplot(p_cfit,xlim,.00001); - title(['PCHIP-interp of ' regexprep(file_name{m},'_','\\_')]); - legend('data','interpolant'); - - subplot(2,1,2); - plot(tube_struct.(t_labels{m})(:,1),tube_struct.(t_labels{m})(:,2),'r*'); - plot(p_cfit,'integral'); - title(['Integral of ' regexprep(file_name{m},'_','\\_')]) - end - end -else -end -assignin('base', 'u_times', u_times); -assignin('base', 'kbq_ml_dkc_mx', kbq_ml_dkc_mx); -assignin('base', 'tube_code_mx', tube_code_mx); - -%% Construct Histograms for extra analysis -if (~population_pf) - make_histograms(u_times, kbq_ml_dkc_mx, tube_code_mx); -end - -%% Fit data -pos_wksp_file=ls(fullfile(home_dir,'PET','Sum','*_workspace.mat')); -if ~isempty(pos_wksp_file) && exist(fullfile(home_dir,'PET','Frame1.lst.hdr'),'file') - wkspfil=fullfile(home_dir,'PET','Sum',deblank(pos_wksp_file(1,:))); - load(deblank(wkspfil(1,:)),'Series_Duration'); - if length(Series_Duration)>tube_struct.(t_labels{m})(end,1); - extrap_time=length(Series_Duration); - else - extrap_time=tube_struct.(t_labels{m})(end,1); - end -else - extrap_time=tube_struct.(t_labels{m})(end,1); -end - -% Fit Parent Fraction according to four different models -[met_fit, GOF_TABLE, fits] = pf_fit('parentfraction.bld','',extrap_time, 4); -assignin('base', 'GOF_TABLE', GOF_TABLE); -assignin('base', 'pf_fits', fits); -set(gcf, 'Position', [50 50 840 630]); -saveas(gcf, 'parentfraction_fits.bmp'); - -% Fit the raw data -update_waitbar(0.7, h, 'Fitting Raw Data (This may take a minute.)'); -raw_bld_file={'plasma_art.bld','plasma_ven.bld','whole-blood_art.bld','whole-blood_ven.bld'}; -for m=1:4 - if exist(fullfile(pwd,raw_bld_file{m}),'file'), - test=dlmread(raw_bld_file{m},'\t',1,0); - if test(1,1)<200 %test(1,1)<60 %only fit data that looks like the peak was captured. - if exist(fullfile(pwd,'parentfraction.bld'),'file') - % feng_fit3 fits data according to Feng & L3Exp models - [bld_fit, models_gof, models]=feng_fit3(raw_bld_file{m},extrap_time,met_fit); - else - [bld_fit, models_gof, models]=feng_fit3(raw_bld_file{m},extrap_time); - end - if exist(fullfile(pwd,'parentfraction_fit.bld'),'file') %%metabolite correction - corr_bld_crv=met_fit.*bld_fit; - end - end - end -end -assignin('base', 'models_gof', models_gof); - -update_waitbar(0.8, h, 'Generating Additional Figures'); -% Use Parent Fraction Fit to apply Met correction (if metcor was selected) -disp 'FENG & LIN3EXP FITS via Lin2Exp Met Correction'; -if exist(fullfile(pwd,'parentfraction_lin2exp.bld'),'file') - bld_fls=cellstr(ls('*.bld')); - bld_fls(cellfun(@(x) ~isempty(x),regexp(bld_fls,'parentfraction','match')))=[]; - bld_fls(cellfun(@(x) ~isempty(x),regexp(bld_fls,'plasma_art.bld','match')))=[]; - bld_fls(cellfun(@(x) ~isempty(x),regexp(bld_fls,'plasma_ven.bld','match')))=[]; - bld_fls(cellfun(@(x) ~isempty(x),regexp(bld_fls,'whole-blood_art.bld','match')))=[]; - bld_fls(cellfun(@(x) ~isempty(x),regexp(bld_fls,'whole-blood_ven.bld','match')))=[]; - bld_fls(cellfun(@(x) ~isempty(x),regexp(bld_fls,'_metcor','match')))=[]; - - PFfit_fls=cellstr(ls('*.bld')); - PFfit_fls(cellfun(@(x) ~isempty(x),regexp(PFfit_fls,'Feng','match')))=[]; - PFfit_fls(cellfun(@(x) ~isempty(x),regexp(PFfit_fls,'L3','match')))=[]; - PFfit_fls(cellfun(@(x) ~isempty(x),regexp(PFfit_fls,'parentfraction.bld','match')))=[]; - PFfit_fls(cellfun(@(x) ~isempty(x),regexp(PFfit_fls,'plasma_art.bld','match')))=[]; - PFfit_fls(cellfun(@(x) ~isempty(x),regexp(PFfit_fls,'plasma_ven.bld','match')))=[]; - PFfit_fls(cellfun(@(x) ~isempty(x),regexp(PFfit_fls,'whole-blood_art.bld','match')))=[]; - PFfit_fls(cellfun(@(x) ~isempty(x),regexp(PFfit_fls,'whole-blood_ven.bld','match')))=[]; - PFfit_fls(cellfun(@(x) ~isempty(x),regexp(PFfit_fls,'linear','match')))=[]; - - iteration = 1; - if ~isempty(bld_fls) - for m=1:size(bld_fls,1) - for i=1:size(PFfit_fls,1) - switch PFfit_fls{i} - case {'parentfraction_hill.bld'} - fname = 'hill'; - case {'parentfraction_exp.bld'} - fname = 'exp'; - case {'parentfraction_power.bld'} - fname = 'power'; - case {'parentfraction_lin2exp.bld'} - fname = 'lin2exp'; - otherwise - error ('Fit model not known to this program.') - end - %iteration = (((m.*size(PFfit_fls,1)) + i) - size(PFfit_fls,1)); - type = metcor_bld2(bld_fls{m},PFfit_fls{i},fname,batoffset,iteration,(size(bld_fls,1)),(size(PFfit_fls,1))); - %metcor_bld3(bld_fls{m},PFfit_fls{i},fname,batoffset); - %metcor_bld4(bld_fls{m},PFfit_fls{i},fname,batoffset);end - if iteration < size(PFfit_fls,1) - iteration = iteration + 1; - else - iteration = 1; - saveas(gcf, fullfile(pwd,strcat('metcor_graphs_',type,'.bmp'))); - end - end - end - end -end - -%% Write header -fid1=fopen(fullfile(csv_path,'blood.hdr'), 'w'); -fprintf(fid1,['TemplateFile:=' regexprep(template_file,'\\','\\\\'), '\r']); -fprintf(fid1,['CSVPath:=' regexprep(csv_path,'\\','\\\\'), '\r']); -fprintf(fid1,['Injection Start:=' num2str(InjectionStart), '\r']); -fprintf(fid1,['BAT Offset:=' num2str(batoffset), '\r']); -fprintf(fid1,['Tracer:=' num2str(tracer), '\r']); -fprintf(fid1,['MetCorrection:=' num2str(met_correct), '\r']); -fprintf(fid1,['Duration:=' num2str(extrap_time), '\r']); -fprintf(fid1,['Reconstruction Date:=' datestr(now), '\r']); -fclose(fid1); - -update_waitbar(0.9, h, 'Creating the Final Report'); -%% Create Final Report (final_report.pdf) -% Assign necessary variables to workspace -assignin('base', 'participant', template_file); -assignin('base', 'dose_info', dose_info); -assignin('base', 'csv_path', csv_path); -assignin('base', 'batoffset', batoffset); -assignin('base', 'num_deleted', num_deleted); -assignin('base', 'models', models); - -if strcmp(report_format, 'pdf') - publish('final_report.m', 'showCode', false, 'maxHeight', 100, 'maxWidth', 100 ... - , 'imageFormat', 'jpg', 'format', 'pdf', 'outputDir', csv_path); -elseif strcmp(report_format, 'html') - publish('final_report.m', 'showCode', false, 'maxHeight', 100, 'maxWidth', 100 ... - , 'imageFormat', 'jpg', 'format', 'html', 'outputDir', csv_path); -elseif strcmp(report_format, 'both') - publish('final_report.m', 'showCode', false, 'maxHeight', 100, 'maxWidth', 100 ... - , 'imageFormat', 'jpg', 'format', 'html', 'outputDir', csv_path); - publish('final_report.m', 'showCode', false, 'maxHeight', 100, 'maxWidth', 100 ... - , 'imageFormat', 'jpg', 'format', 'pdf', 'outputDir', csv_path); -end - -%% Clean Up -update_waitbar(0.95, h, 'Cleaning Up'); -% Put HTML report, BMP Images, & PNG Images into a separate folder -mkdir 'images_and_html'; -if strcmp(report_format, 'html') || strcmp(report_format, 'both') - copyfile('final_report.html', 'images_and_html'); - delete('final_report.html'); -end -bmps = ls('*.bmp'); -for m=1:size(bmps) - file = regexprep(bmps(m,:), ' ', ''); - %disp( [file, '.']); - copyfile(file, 'images_and_html'); - delete(file); -end -pngs = ls('*.png'); -for m=1:size(pngs) - file = regexprep(pngs(m,:), ' ', ''); - %disp( [file, '.']); - copyfile(file, 'images_and_html'); - delete(file); -end -addpath('images_and_html'); - -% Put BLD files into a separate folder -mkdir 'bld_files'; -blds = ls('*.bld'); -for m=1:size(blds) - file = regexprep(blds(m,:), ' ', ''); - %disp( [file, '.']); - copyfile(file, 'bld_files'); - delete(file); -end -addpath('bld_files'); - -update_waitbar(1, h, 'Analysis Complete!'); - -% Clear Workspace Variables -clear; - -disp ' *****ANALYSIS COMPLETE***** '; -close all; -end - - - -%% ======================================================================= -function [tube_struct,t_labels]=sort_tube_codes(tube_code_mx,col_time_mx,kbq_ml_dkc_mx) -% count_time_mx -% count_dur_mx -% t_counts_mx -% dev_cpm_mx -% t_counts_error_mx -% counts_dkc_mx -%% First, sort tubes by label then time -t_labels=unique(tube_code_mx(:)); -t_labels(cellfun(@(x) any(strcmpi(x,{'B','E','K'})),t_labels))=[]; %remove B and E -tube_struct=struct; -for m=1:size(t_labels,1) - lab_pos=find(strcmpi(tube_code_mx,t_labels{m})); - [time,lab_pos_ind]=sort(col_time_mx(lab_pos)); - ordered_data=kbq_ml_dkc_mx(lab_pos(lab_pos_ind)); - tube_struct.(t_labels{m})=[time, ordered_data]; - tube_struct.(t_labels{m})(tube_struct.(t_labels{m})(:,1)<0,:)=[]; %remove negative times - tube_struct.(t_labels{m})(tube_struct.(t_labels{m})(:,2)<0,:)=[]; %remove negative numbers -end -end - - - -%% ======================================================================= -function exit_gamma_bomb = check_overwrite(csv_path, population_pf) -% Input: the directory of csv_files -% Output: Checks to make sure a final report has not already been generated -% for this data. If a final report is found: -% - Delete existing final report -% - Delete any exisiting image files/graphs that were generated before -% - Delete any existing bld files that were generated beforeif (exist(fullfile(csv_path,'final_report.html'),'file')) - images=ls(fullfile(csv_path,'*.bmp')); - disp(images); - [num_imgs,~] = size(images); - - overwrite_warning = 'Empty'; - if ( exist([csv_path '\final_report.pdf'], 'file') || exist([csv_path '\final_report.html'], 'file') || (num_imgs>0) ||... - exist([csv_path '\bld_files'], 'dir') || exist([csv_path '\images_and_html'], 'dir') ) - overwrite_warning = questdlg(... - 'Final Report or images already generated for this participant. Overwrite? (All existing final reports, bld files, and bitmap images will be deleted.)'... - ,'Overwrite Warning','Yes','No','No'); - end - if strcmp(overwrite_warning, 'No') - disp 'Program Terminated.'; - exit_gamma_bomb = 1; - elseif strcmp(overwrite_warning, 'Empty') - exit_gamma_bomb = 0; - return; - else - if exist(fullfile(csv_path,'final_report.html'),'file') - delete(fullfile(csv_path,'final_report.html'),'file'); - end - if exist(fullfile(csv_path,'final_report.pdf'),'file') - delete(fullfile(csv_path,'final_report.pdf'),'file'); - end - if exist([csv_path '\bld_files'], 'dir') - rmdir('bld_files', 's'); - end - if exist([csv_path '\images_and_html'], 'dir') - rmdir('images_and_html', 's'); - end - images=ls(fullfile(csv_path,'*.bmp')); - pngs = ls(fullfile(csv_path,'*.png')); - images = cat(1, images, pngs); - [num_img,~] = size(images); - disp(num_img); - if (num_img > 0) - for m=1:num_img - delete(strcat(csv_path,'\',images(m,:))); - end - end - - % Don't delete bld files if you're using a population parent - % fraction! - if (~population_pf) - builds=ls(fullfile(csv_path,'*.bld')); - - [num_bld,~] = size(builds); - disp(num_bld); - if (num_bld > 0) - for m=1:num_bld - delete(strcat(csv_path,'\',builds(m,:))); - end - end - end - - - clear images builds - exit_gamma_bomb = 0; - %delete(images); - %delete(builds); - end -end - - -%% ======================================================================== -function update_waitbar(val, h, message) -% Updates the waitbar to show program's progress -if val == 0 - h = waitbar(0, [message, ' - ', num2str(val*100), '%']); -elseif val > 1 - return; -elseif val == 1 - close(h); -else - waitbar(val, h, [message, ' - ', num2str(val*100), '%']); -end -end - diff --git a/gamma_bomb_2/code/gamma_bomb_gui.fig b/gamma_bomb_2/code/gamma_bomb_gui.fig deleted file mode 100644 index 29ac3f0..0000000 Binary files a/gamma_bomb_2/code/gamma_bomb_gui.fig and /dev/null differ diff --git a/gamma_bomb_2/code/gamma_bomb_gui.m b/gamma_bomb_2/code/gamma_bomb_gui.m deleted file mode 100644 index ce18b27..0000000 --- a/gamma_bomb_2/code/gamma_bomb_gui.m +++ /dev/null @@ -1,373 +0,0 @@ -function varargout = gamma_bomb_gui(varargin) -% GAMMA_BOMB_GUI MATLAB code for gamma_bomb_gui.fig -% GAMMA_BOMB_GUI, by itself, creates a new GAMMA_BOMB_GUI or raises the existing -% singleton*. -% -% H = GAMMA_BOMB_GUI returns the handle to a new GAMMA_BOMB_GUI or the handle to -% the existing singleton*. -% -% GAMMA_BOMB_GUI('CALLBACK',hObject,eventData,handles,...) calls the local -% function named CALLBACK in GAMMA_BOMB_GUI.M with the given input arguments. -% -% GAMMA_BOMB_GUI('Property','Value',...) creates a new GAMMA_BOMB_GUI or raises the -% existing singleton*. Starting from the left, property value pairs are -% applied to the GUI before gamma_bomb_gui_OpeningFcn gets called. An -% unrecognized property name or invalid value makes property application -% stop. All inputs are passed to gamma_bomb_gui_OpeningFcn via varargin. -% -% *See GUI Options on GUIDE's Tools menu. Choose "GUI allows only one -% instance to run (singleton)". -% -% See also: GUIDE, GUIDATA, GUIHANDLES - -% Edit the above text to modify the response to help gamma_bomb_gui - -% Last Modified by GUIDE v2.5 03-Jun-2015 15:19:03 - -% Begin initialization code - DO NOT EDIT -gui_Singleton = 1; -gui_State = struct('gui_Name', mfilename, ... - 'gui_Singleton', gui_Singleton, ... - 'gui_OpeningFcn', @gamma_bomb_gui_OpeningFcn, ... - 'gui_OutputFcn', @gamma_bomb_gui_OutputFcn, ... - 'gui_LayoutFcn', [] , ... - 'gui_Callback', []); -if nargin && ischar(varargin{1}) - gui_State.gui_Callback = str2func(varargin{1}); -end - -if nargout - [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:}); -else - gui_mainfcn(gui_State, varargin{:}); -end -% End initialization code - DO NOT EDIT - - -% --- Executes just before gamma_bomb_gui is made visible. -function gamma_bomb_gui_OpeningFcn(hObject, eventdata, handles, varargin) -% This function has no output args, see OutputFcn. -% hObject handle to figure -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) -% varargin command line arguments to gamma_bomb_gui (see VARARGIN) - -% Choose default command line output for gamma_bomb_gui -handles.output = hObject; - -% Update handles structure -guidata(hObject, handles); - -% UIWAIT makes gamma_bomb_gui wait for user response (see UIRESUME) -% uiwait(handles.gamma_gui); - - -% --- Outputs from this function are returned to the command line. -function varargout = gamma_bomb_gui_OutputFcn(hObject, eventdata, handles) -% varargout cell array for returning output args (see VARARGOUT); -% hObject handle to figure -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Get default command line output from handles structure -varargout{1} = handles.output; - - - -function xls_ffile_Callback(hObject, eventdata, handles) -% hObject handle to xls_ffile (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hints: get(hObject,'String') returns contents of xls_ffile as text -% str2double(get(hObject,'String')) returns contents of xls_ffile as a double - - -% --- Executes during object creation, after setting all properties. -function xls_ffile_CreateFcn(hObject, eventdata, handles) -% hObject handle to xls_ffile (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: edit controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); -end - - -% --- Executes on button press in get_xls_ffile. -function get_xls_ffile_Callback(hObject, eventdata, handles) -% hObject handle to get_xls_ffile (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) -[filename,filepath,~]=uigetfile(fullfile(pwd,'*.xls;*xlsx'),'Pick Template File'); -cd(filepath) -set(handles.xls_ffile,'String',fullfile(filepath,filename)) -guidata(hObject, handles); - - -function csv_dir_Callback(hObject, eventdata, handles) -% hObject handle to csv_dir (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hints: get(hObject,'String') returns contents of csv_dir as text -% str2double(get(hObject,'String')) returns contents of csv_dir as a double - - -% --- Executes during object creation, after setting all properties. -function csv_dir_CreateFcn(hObject, eventdata, handles) -% hObject handle to csv_dir (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: edit controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); -end - - -% --- Executes on button press in get_csv_dir. -function get_csv_dir_Callback(hObject, eventdata, handles) -% hObject handle to get_csv_dir (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) -set(handles.csv_dir,'String',uigetdir(pwd)) -guidata(hObject, handles); - -% --- Executes on selection change in tracer. -function tracer_Callback(hObject, eventdata, handles) -% hObject handle to tracer (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hints: contents = cellstr(get(hObject,'String')) returns tracer contents as cell array -% contents{get(hObject,'Value')} returns selected item from tracer - - -% --- Executes during object creation, after setting all properties. -function tracer_CreateFcn(hObject, eventdata, handles) -% hObject handle to tracer (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: popupmenu controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); -end - - -% --- Executes on button press in allons_y. -function allons_y_Callback(hObject, eventdata, handles) -% hObject handle to allons_y (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) -assignin('base','ligand',get(handles.ligand,'String')); -cd(get(handles.csv_dir,'String')); - -% Quickly check data to make sure all the times make sense -dlg={'***Template***'}; -template_TOI=datestr(xlsread(get(handles.xls_ffile,'String'),4),'HH:MM:SS'); -dlg(end+1)={['Template TOI = ' template_TOI]}; -col_mx=xlsread(get(handles.xls_ffile,'String'),2); -first_sample_t=datestr(col_mx(2,1),'HH:MM:SS'); -dlg(end+1)={['First Draw Time = ' first_sample_t]}; - -% First csv file -dlg(end+1)={'***CSV***'}; -csv_files=ls(get(handles.csv_dir,'String')); -fid_in = fopen(fullfile(get(handles.csv_dir,'String'),csv_files(3,:)),'r'); -csv_dat_str = '%d%s%s%d%d%d%d%d%f%s%f%f%f%s'; -csv_dat= textscan(fid_in,csv_dat_str,'delimiter',',','HeaderLines',1,'EndOfLine','\n'); -csv_read_t=datestr(csv_dat{3}{1},'HH:MM:SS'); -dlg(end+1)={['First Read Time = ' csv_read_t]}; - -% Bay 6 Dose info and Frame1.hdr -hd=home_dir(get(handles.xls_ffile,'String')); -if exist(fullfile(hd,'PET','Dose_info.xls'),'file') - dlg(end+1)={'***Dose Info***'}; - di=xlsread(fullfile(hd,'PET','Dose_info.xls'),1); - dose_calib_t=datestr(di(3,7),'HH:MM:SS'); - dlg(end+1)={['Dose Calib. Time = ' dose_calib_t]}; - Series_t=datestr(di(3,1),'HH:MM:SS'); - dlg(end+1)={['Series Time = ' Series_t]}; - Acq_t=datestr(di(3,1),'HH:MM:SS'); - dlg(end+1)={['Acq Time = ' Acq_t]}; -end -dlg(end+1)={'***Frame 1***'}; -if exist(fullfile(hd,'PET','Frame1.lst.hdr'),'file') - pd=read_i_hdr(fullfile(hd,'PET','Frame1.lst.hdr')); - PET_start_t=pd.StudyTimeHhMmSs; -elseif exist(fullfile(hd,'PET','Frame1.hdr'),'file') - pd=read_i_hdr(fullfile(hd,'PET','Frame1.hdr')); - PET_start_t=pd.StudyTimeHhMmSs; -else - PET_start_t=''; -end -dlg(end+1)={['Frame1.lst.hdr Start Time = ' PET_start_t]}; -dlg(end+1)={''}; -dlg(end+1)={'Do these times make sense?'}; -button=questdlg(dlg); - -% Was a population curve generated for the parent fraction? -population_bool = false; -if (strcmp(get(handles.pop_gen_button, 'String'), 'Parent Fraction Ready!')) - population_bool = true; -end - -% Begin Data Analysis if times make sense -switch button - case 'Yes' -contents = cellstr(get(handles.tracer,'String')); -format = cellstr(get(handles.report_format,'String')); - display(get(handles.xls_ffile,'String')) - display(get(handles.csv_dir,'String')) - display(contents{get(handles.tracer,'Value')}) - display(str2num(get(handles.bkgnd,'String'))) -if isempty(get(handles.bkgnd,'String')) -gamma_bomb(get(handles.xls_ffile,'String'),get(handles.csv_dir,'String'),... - contents{get(handles.tracer,'Value')},get(handles.met_cor,'Value'),... - get(handles.dose_info_ffile, 'String'), format{get(handles.report_format, 'Value')},population_bool) -else - display('Static Background Entered...') -gamma_bomb(get(handles.xls_ffile,'String'),get(handles.csv_dir,'String'),... - contents{get(handles.tracer,'Value')},get(handles.met_cor,'Value'),... - get(handles.dose_info_ffile, 'String'), format{get(handles.report_format, 'Value')},population_bool,str2num(get(handles.bkgnd,'String'))) -end -otherwise -end - -function bkgnd_Callback(hObject, eventdata, handles) -% hObject handle to bkgnd (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hints: get(hObject,'String') returns contents of bkgnd as text -% str2double(get(hObject,'String')) returns contents of bkgnd as a double - - -% --- Executes during object creation, after setting all properties. -function bkgnd_CreateFcn(hObject, eventdata, handles) -% hObject handle to bkgnd (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: edit controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); -end - - -% --- Executes on button press in met_cor. -function met_cor_Callback(hObject, eventdata, handles) -% hObject handle to met_cor (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hint: get(hObject,'Value') returns toggle state of met_cor - - -% --- Executes on button press in adj_csv. -function adj_csv_Callback(~, eventdata, handles) -% hObject handle to adj_csv (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) -time_inc_GUI - - -% --- Executes on button press in dose_info. -function dose_info_Callback(hObject, eventdata, handles) -% hObject handle to dose_info (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) -[filename,filepath,~]=uigetfile(fullfile(pwd,'*.xls;*xlsx'),'Pick Dose Info File'); -cd(filepath) -set(handles.dose_info_ffile,'String',fullfile(filepath,filename)) -guidata(hObject, handles); - - - -function dose_info_ffile_Callback(hObject, eventdata, handles) -% hObject handle to dose_info_ffile (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hints: get(hObject,'String') returns contents of dose_info_ffile as text -% str2double(get(hObject,'String')) returns contents of dose_info_ffile as a double - - -% --- Executes during object creation, after setting all properties. -function dose_info_ffile_CreateFcn(hObject, eventdata, handles) -% hObject handle to dose_info_ffile (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: edit controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); -end - - -% --- Executes on selection change in report_format. -function report_format_Callback(hObject, eventdata, handles) -% hObject handle to report_format (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hints: contents = cellstr(get(hObject,'String')) returns report_format contents as cell array -% contents{get(hObject,'Value')} returns selected item from report_format - - -% --- Executes during object creation, after setting all properties. -function report_format_CreateFcn(hObject, eventdata, handles) -% hObject handle to report_format (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: popupmenu controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); -end - - -% --- Executes on button press in pop_gen_button. -function pop_gen_button_Callback(hObject, eventdata, handles) -% hObject handle to pop_gen_button (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) -if (strcmp(get(handles.csv_dir,'String'),'Edit Text')) - warndlg('Please choose a CSV dir before continuing.', 'Must Select CSV Directory'); -else - pop_gen_gui(); - %set(handles.pop_gen_button, 'Enable', 'off', 'String', 'Parent Fraction Ready!'); -end - - - -function ligand_Callback(hObject, eventdata, handles) -% hObject handle to ligand (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hints: get(hObject,'String') returns contents of ligand as text -% str2double(get(hObject,'String')) returns contents of ligand as a double - - -% --- Executes during object creation, after setting all properties. -function ligand_CreateFcn(hObject, eventdata, handles) -% hObject handle to ligand (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: edit controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); - -end diff --git a/gamma_bomb_2/code/generate_head_curve.m b/gamma_bomb_2/code/generate_head_curve.m deleted file mode 100644 index b662a5e..0000000 --- a/gamma_bomb_2/code/generate_head_curve.m +++ /dev/null @@ -1,56 +0,0 @@ -function [head_curve,head_curve_trues]=generate_head_curve(file,Series_Duration) - -acq_timing_error=0; - -head_curve_full = csvread(file, 1, 0); -if head_curve_full(1,1)~=1000 - warndlg(sprintf('Head curve file was found to start at a time~=1 sec. \n I''m going to attempt to fix this.')); - time_offset=head_curve_full(1,1); - head_curve_full(:,1)=head_curve_full(:,1)-head_curve_full(1,1)+1000; - acq_timing_error=1; -end - -if nargin<2, Series_Duration(1,1:size(head_curve_full,1))=1;end - -head_curve_reduced=head_curve_full(:,1:2); -head_curve_reduced(:,1)=floor(head_curve_reduced(:,1)/1000); - -% compute trues as ERS_prompts - ERS_delays -head_curve_true_reduced=[head_curve_full(:,1) head_curve_full(:,10)-head_curve_full(:,11)]; -head_curve_true_reduced(head_curve_true_reduced<0)=NaN; - -head_curve_true_reduced(:,1)=floor(head_curve_true_reduced(:,1)/1000); -head_curve=Series_Duration'; -head_curve_trues=Series_Duration'; - -%%align Head Curve to Number Line -for m=1:size(head_curve_reduced,1); - if head_curve_reduced(m,1)=ts & t=ts & t=ts+tp; - -aif(pre)=0; -aif(inc)=A1*(t(inc)-ts); -aif(post)=(A1*(tp)+A2+A3-C)/3*exp(-L1*(t(post)-(ts+tp)))+(A1*(tp)-A2-C)/3*exp(-L2*(t(post)-(ts+tp)))+(A1*(tp)-A3-C)/3*exp(-L3*(t(post)-(ts+tp)))+C; - - -end \ No newline at end of file diff --git a/gamma_bomb_2/code/linear_2exp_model.m b/gamma_bomb_2/code/linear_2exp_model.m deleted file mode 100644 index 3cbc980..0000000 --- a/gamma_bomb_2/code/linear_2exp_model.m +++ /dev/null @@ -1,15 +0,0 @@ -function aif=linear_2exp_model(A1,A2,A3,tau,lambda1,lambda2,t) -%% -%Y(t)=A1*t+A2 0=tau; - -aif(pre)=A1*t(pre)+A2; -aif(post)=(A1*tau+A2+A3)/2*exp(-lambda1*(t(post)-tau))+(A1*tau+A2-A3)/2*exp(-lambda2*(t(post)-tau)) ; - -end \ No newline at end of file diff --git a/gamma_bomb_2/code/make_histograms.m b/gamma_bomb_2/code/make_histograms.m deleted file mode 100644 index 785b9f3..0000000 --- a/gamma_bomb_2/code/make_histograms.m +++ /dev/null @@ -1,92 +0,0 @@ -function make_histograms(times, kbq_ml_dkc_mx, t_codes) -% Creates a 9-bar histogram for each time point in the parent fraction data -% times = Titles (each histogram represents one time) -% kbq_ml_dkc_mx = y-data -% t_codes = matrix of data types (same dimensions as kbq_ml_dkc_mx) ('M = -% Metabolite', 'P = Parent', etc...) - -% Get dimensions of input matrices -[rows, cols] = size(t_codes); - -% Figure out how many subplots (num_graphs) will be needed and which row of -% the kbq_ml_dkc matrix will be graphed as the first subplot (first) -first = rows; -num_graphs = 0; -for r=1:rows - will_graph = false; - for c=1:cols - if strcmp(t_codes(r,c),'M') || strcmp(t_codes(r,c),'P') - if r < first - first = r-1; - end - will_graph = true; - end - end - if will_graph - num_graphs = num_graphs + 1; - end -end - -% Plot the data -for iter=1:2 - figure; - hold on; - for r=1:rows - sub = false; - for c=1:cols - if strcmp(t_codes(r,c),'M') || strcmp(t_codes(r,c),'P') - sub = true; - end - end - if sub - subplot(1,num_graphs,r-first); - met = zeros(1,0); - parent = zeros(1,0); - num_m_bars = 0; - num_p_bars = 0; - for c=1:cols - if strcmp(t_codes(r,c),'M') - met(end+1)=kbq_ml_dkc_mx(r,c); - num_m_bars = num_m_bars + 1; - elseif strcmp(t_codes(r,c), 'P') - parent(end+1)=kbq_ml_dkc_mx(r,c); - num_p_bars = num_p_bars+1; - end - end - % Put the parent & met data into different columns to allow for - % different colors - data = zeros(num_m_bars+num_p_bars,2); - data(1:num_m_bars,1) = met'; - data(num_m_bars+1:end,2) = parent'; - h = bar(data,3); - set(h(1),'FaceColor','r'); - set(h(2),'FaceColor',[.2, 0, 0.8]); - if iter == 1 - ylim([0 1]); - xlabel('Tube #'); - ylabel('kBq/mL'); - else - ylim('auto'); - xlabel('Tube #'); - ylabel('kBq/mL'); - end - if (r-first == 1) - legend('Metabolite', 'Parent', 'Location', 'NorthWestOutside'); - else - legend 'off'; - end - title_str = [num2str(round(times(r-first)./60)) ' Min.']; - title(title_str, 'FontSize', 10, 'FontWeight', 'bold'); - - set(gca, 'FontSize', 7); - end - end - -% Save Histograms as a BMP image -set(gcf, 'Position', [519 711 1250 270]); -disp(iter) -filename = ['histograms_' num2str(iter)]; -saveas(gcf, filename, 'bmp'); -end - -end \ No newline at end of file diff --git a/gamma_bomb_2/code/metcor_bld2.m b/gamma_bomb_2/code/metcor_bld2.m deleted file mode 100644 index e2481b7..0000000 --- a/gamma_bomb_2/code/metcor_bld2.m +++ /dev/null @@ -1,99 +0,0 @@ -%% MODIFIED by Tom Morin (05/01/2015) from script by Daniel Chonde - -function type = metcor_bld2(blood_ffile,met_ffile,fit_source,batoffset,iteration,figrows,figcols) - % This creates the various metabolite correction bld files & BAT - % correction bld files. - % INPUT: - % - blood_ffile = File containing fitted TAC blood data (i.e plasma_art_Feng.bld) - % - met_ffile = File containing Parent Fraction data (i.e. parent_fraction_hill.bld) - % - fit_source = name of the Parent Fraction fit (i.e. hill, exp, etc.) - % - batoffset = Bolus Arrival Time offset in seconds (calculated previously) - % - Iteration, figrows, & figcols, were all used to organize the - % aesthetics of the subplots/figures. These are no longer needed in - % the final report, so they are no longer produced. - % OUTPUT: - % - bld files for the met correction & BAT correction of each model fit - % - type = A name describing the model fit & its met correciton (i.e. Feng_Hill, L3Exp_Power, etc.) - - - % Read in data & create necessary matrices - raw_data=dlmread('plasma_art.bld','\t',1,0); - plasma_art=dlmread(blood_ffile,'\t',1,0); - parentfraction=dlmread(met_ffile,'\t',1,0); - - %This method for finding sigfigs doesn't work when the last value is a - %whole number -% teststr=num2str(plasma_art(find(plasma_art(:,2)>0,1,'last'),2)); -% splits=regexp(teststr,'\.','split'); -% sf=size(splits{2},2); - fid=fopen(blood_ffile,'r'); - fgetl(fid); - teststr=fgetl(fid); - splits=regexp(teststr,'\.','split'); - sf=length(splits{end}); - fclose(fid); - - if length(plasma_art(:,1))==length(parentfraction(:,1)) - plasma_art_metcor(:,1)=plasma_art(:,1); - plasma_art_metcor(:,2)=plasma_art(:,2).*parentfraction(:,2); - else - error('Dimensions do not match'); - end - - % Set file name - name = regexprep(blood_ffile, '.bld', '_metcor'); - name = strcat(name, '_', fit_source); - type = regexprep(name, 'plasma_art_', ''); - type = regexprep(type, strcat('_metcor_', fit_source), ''); - -% * * The Plots generated below are no longer used in the final report. * * -% ========================================================================= -% % Plot some shit -% if iteration == 1 -% figure -% end -% subplot(figrows.*2,figcols./2,((iteration.*2)-1));plot(plasma_art(:,1),plasma_art(:,2));hold on; -% plot(plasma_art_metcor(:,1),plasma_art_metcor(:,2),'r'); -% %plot(raw_data(:,1),raw_data(:,2),'xk'); -% set(gca,'XScale','log', 'FontSize', 7); -% -% ylabel('kBq/ml', 'FontSize', 7); -% if iteration == 1 -% legend({'uncorrected','met corrected'}, 'Location', 'NorthEast', 'FontSize', 6, 'Orientation', 'vertical') -% title(strcat(type, ' Fits'), 'FontSize', 12, 'FontWeight', 'bold') -% else -% legend('off'); -% end -% subplot(figrows.*2,figcols./2,(iteration.*2));plot(parentfraction(:,1),parentfraction(:,2), 'm') -% title(strcat(upper(fit_source), ' Fit Parent-Fraction'), 'FontSize', 9, 'FontWeight', 'bold'); -% ylabel('Fraction', 'FontSize', 7); -% xlabel('Time (s)', 'FontSize', 7); -% set(gca, 'FontSize', 7); -% %hold off; -% -% set(gcf, 'Position', [50 50 750 800]); -% ========================================================================= -% * * The Plots generated above are no longer used in the final report. * * - - % Save JPG image of graph - % saveas(gcf,fullfile(pwd,strcat(name, '.bmp'))); - % set(gca,'XScale','log'); - - % Save BLD file with fit data - names=textread(blood_ffile,'%s',1,'delimiter','\n'); - hdrnames=regexp(names,'\t','split'); - print_list(hdrnames{:},plasma_art_metcor,strcat(name, '.bld'),sf); - - % Apply BAT Offset & Save BLD file - A=zeros(batoffset,2); - [len,~]=size(A); - for m=1:len-1 - A(m,1) = m; - end - for m=1:length(plasma_art_metcor) - plasma_art_metcor(m,1) = plasma_art_metcor(m,1) + batoffset; - end - bat_cor = cat(1, A, plasma_art_metcor); - print_list(hdrnames{:},bat_cor,strcat(name, '_BATcor.bld'),sf); -end - diff --git a/gamma_bomb_2/code/model_aif.m b/gamma_bomb_2/code/model_aif.m deleted file mode 100644 index dcd02d1..0000000 --- a/gamma_bomb_2/code/model_aif.m +++ /dev/null @@ -1,20 +0,0 @@ -function mu=model_aif(a,b,c,b1,b2,t) -%%Model of AIF bolus for extraction of bolus arrival time (BAT) -% MU=c, ta+b -%% at some point I need to come back to this and change the equation so b is the amount of time where the first line is used -%since b must always be greater than a. This can be enforced by writing -%the equation such that a<=t<=a+b, or essentially replacing b by (b+a) - -mu=zeros(size(t)); - -pre= a>=t; -low_side= (a std_cutoff*std(met); - outliers = excludedata(time,met,'indices',I); - options.Exclude=outliers; - [fit2,gof2] = fit(time,met,ffun,options); - - % Plot the data - subplot(2,2,iter); - if sum(outliers) > 0 - plot(fit2,'r-', time, met, '.b',outliers,'r*'); - else - plot(fit2, 'r-', time, met, '.b'); - end - hold on; - %plot(fit1, 'b:'); % This is the original fit, before outliers deleted - - % Set up Legend - if iter == 1 && sum(outliers) > 0 - legend('Data', 'Oultiers', 'Fit', 'Location', 'NorthEast' ); - elseif iter == 1 - legend('Data', 'Fit', 'Location', 'NorthEast'); - else - legend 'off'; - end - - % Label axes - xlabel('Time (s)', 'FontSize', 7); - ylabel('Parent Fraction', 'FontSize', 7'); - set(gca, 'FontSize', 7); - title([upper(fit_type) ' Fit of Parent Fraction'], 'FontSize', 10, 'FontWeight', 'bold'); - - % Denote Goodness of Fit Data - GOF_TABLE(iter,1) = gof2.rsquare; - GOF_TABLE(iter,2) = gof2.adjrsquare; - GOF_TABLE(iter,3) = gof2.dfe; - - % Determine Best Fit & Pass it back for met correction - if iter == 1 - best_fit = fit2; - elseif abs(GOF_TABLE(iter,1)-1) < abs(GOF_TABLE(iter-1,1)-1) - best_fit = fit2; - end - - % Save the fitted data to a bld file - fmet=fit2(time(1):extrap_time); - [metpath,metname,metext]=fileparts(met_ffile); - print_list(hdrnames{:},[[time(1):extrap_time]', fmet],fullfile(metpath,strcat(metname,'_',fit_type,metext)),sf); -end -end \ No newline at end of file diff --git a/gamma_bomb_2/code/pop_gen_gui.fig b/gamma_bomb_2/code/pop_gen_gui.fig deleted file mode 100644 index c6494d6..0000000 Binary files a/gamma_bomb_2/code/pop_gen_gui.fig and /dev/null differ diff --git a/gamma_bomb_2/code/pop_gen_gui.m b/gamma_bomb_2/code/pop_gen_gui.m deleted file mode 100644 index 66b60fd..0000000 --- a/gamma_bomb_2/code/pop_gen_gui.m +++ /dev/null @@ -1,341 +0,0 @@ -function varargout = pop_gen_gui(varargin) -% POP_GEN_GUI MATLAB code for pop_gen_gui.fig -% POP_GEN_GUI, by itself, creates a new POP_GEN_GUI or raises the existing -% singleton*. -% -% H = POP_GEN_GUI returns the handle to a new POP_GEN_GUI or the handle to -% the existing singleton*. -% -% POP_GEN_GUI('CALLBACK',hObject,eventData,handles,...) calls the local -% function named CALLBACK in POP_GEN_GUI.M with the given input arguments. -% -% POP_GEN_GUI('Property','Value',...) creates a new POP_GEN_GUI or raises the -% existing singleton*. Starting from the left, property value pairs are -% applied to the GUI before pop_gen_gui_OpeningFcn gets called. An -% unrecognized property name or invalid value makes property application -% stop. All inputs are passed to pop_gen_gui_OpeningFcn via varargin. -% -% *See GUI Options on GUIDE's Tools menu. Choose "GUI allows only one -% instance to run (singleton)". -% -% See also: GUIDE, GUIDATA, GUIHANDLES - -% Edit the above text to modify the response to help pop_gen_gui - -% Last Modified by GUIDE v2.5 29-May-2015 14:33:13 - -% Begin initialization code - DO NOT EDIT -gui_Singleton = 1; -gui_State = struct('gui_Name', mfilename, ... - 'gui_Singleton', gui_Singleton, ... - 'gui_OpeningFcn', @pop_gen_gui_OpeningFcn, ... - 'gui_OutputFcn', @pop_gen_gui_OutputFcn, ... - 'gui_LayoutFcn', [] , ... - 'gui_Callback', []); -if nargin && ischar(varargin{1}) - gui_State.gui_Callback = str2func(varargin{1}); -end - -if nargout - [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:}); -else - gui_mainfcn(gui_State, varargin{:}); -end -% End initialization code - DO NOT EDIT - - -% --- Executes just before pop_gen_gui is made visible. -function pop_gen_gui_OpeningFcn(hObject, eventdata, handles, varargin) -% This function has no output args, see OutputFcn. -% hObject handle to figure -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) -% varargin command line arguments to pop_gen_gui (see VARARGIN) - -% Choose default command line output for pop_gen_gui -handles.output = hObject; - -% Update handles structure -guidata(hObject, handles); - -h = findobj('Tag', 'gamma_gui'); -gamma_bomb_data = guidata(h); -directory = get(gamma_bomb_data.csv_dir, 'String'); -set(handles.place,'String', directory); -% UIWAIT makes pop_gen_gui wait for user response (see UIRESUME) -% uiwait(handles.population_generation); - - -% --- Outputs from this function are returned to the command line. -function varargout = pop_gen_gui_OutputFcn(hObject, eventdata, handles) -% varargout cell array for returning output args (see VARARGOUT); -% hObject handle to figure -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Get default command line output from handles structure -varargout{1} = handles.output; - - - -% --- Executes during object creation, after setting all properties. -function listbox1_CreateFcn(hObject, eventdata, handles) -% hObject handle to listbox1 (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: listbox controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); -end - - -% --- Executes on selection change in popupmenu1. -function popupmenu1_Callback(hObject, eventdata, handles) -% hObject handle to popupmenu1 (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hints: contents = cellstr(get(hObject,'String')) returns popupmenu1 contents as cell array -% contents{get(hObject,'Value')} returns selected item from popupmenu1 - - -% --- Executes during object creation, after setting all properties. -function popupmenu1_CreateFcn(hObject, eventdata, handles) -% hObject handle to popupmenu1 (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: popupmenu controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); -end - - -% --- Executes on selection change in list. -function list_Callback(hObject, eventdata, handles) -% hObject handle to list (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hints: contents = cellstr(get(hObject,'String')) returns list contents as cell array -% contents{get(hObject,'Value')} returns selected item from list - - -% --- Executes during object creation, after setting all properties. -function list_CreateFcn(hObject, eventdata, handles) -% hObject handle to list (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: listbox controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); -end - - -% --- Executes on button press in browse. -function browse_Callback(hObject, eventdata, handles) -% hObject handle to browse (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) -[filename,filepath,~]=uigetfile(fullfile(pwd,'*.bld'),'Pick Data File'); -cd(filepath) -set(handles.cur_file,'String',fullfile(filepath,filename)) -guidata(hObject, handles); - - -% --- Executes on button press in add. -function add_Callback(hObject, eventdata, handles) -% hObject handle to add (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Add a selected file to the population list -entries = get(handles.list,'String'); -value = get(handles.list,'Value'); -newEntryName = { get(handles.cur_file,'String') }; - -if ~strcmp(newEntryName,'Filename') - on_list = false; - for m=1:length(entries) - if strcmp(entries(m),newEntryName) - on_list = true; - end - end - if ~on_list - if value > 0 && strcmp(entries(1,:),'List of bld Data Files') - entries = newEntryName; - value = 1; - else - entries = [entries; newEntryName]; - value = value+1; - end - if value > 0 - set(handles.remove, 'Enable', 'on'); - end - set(handles.list,'String',entries,'Value',value); - else - warndlg('This file has already been added to the list.'); - end -end - - -% --- Executes on button press in remove. -function remove_Callback(hObject, eventdata, handles) -% hObject handle to remove (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Remove a selected file from the population list -entries = get(handles.list, 'String'); -value = get(handles.list, 'Value'); - -entries(value) = []; -nentries = length(entries); - -if value > nentries - value = value-1; -end - -set(handles.list,'Value',value,'String',entries) - -if nentries == 0 - set(handles.remove,'Enable','off') -end - - -function cur_file_Callback(hObject, eventdata, handles) -% hObject handle to cur_file (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hints: get(hObject,'String') returns contents of cur_file as text -% str2double(get(hObject,'String')) returns contents of cur_file as a double - - -% --- Executes during object creation, after setting all properties. -function cur_file_CreateFcn(hObject, eventdata, handles) -% hObject handle to cur_file (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: edit controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); -end - - -% --- Executes on button press in analyze. -function analyze_Callback(hObject, eventdata, handles) -% hObject handle to analyze (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Produce Parent Fraction Data from the selected population and continue on -% to analyze the rest of the data (returns to Gamma_Bomb_GUI -complete = false; -if strcmp(get(handles.list,'String'),'List of bld Data Files') - warndlg('Please add at least one bld file to the list of data files.', 'File List Empty'); -else -complete = pop_parent_fraction(get(handles.list,'String'), get(handles.place,'String')); -end -if complete - h = findobj('Tag', 'gamma_gui'); - gamma_bomb_data = guidata(h); - set(gamma_bomb_data.pop_gen_button, 'Enable', 'off', 'String', 'Parent Fraction Ready!'); - close(pop_gen_gui); - return; -end - - - -function place_Callback(hObject, eventdata, handles) -% hObject handle to place (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Hints: get(hObject,'String') returns contents of place as text -% str2double(get(hObject,'String')) returns contents of place as a double - - -% --- Executes during object creation, after setting all properties. -function place_CreateFcn(hObject, eventdata, handles) -% hObject handle to place (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles empty - handles not created until after all CreateFcns called - -% Hint: edit controls usually have a white background on Windows. -% See ISPC and COMPUTER. -if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) - set(hObject,'BackgroundColor','white'); -end - - -% --- Executes on button press in open_pop. -function open_pop_Callback(hObject, eventdata, handles) -% hObject handle to open_pop (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) - -% Read in the imported population -entries = get(handles.list, 'String'); -value = get(handles.list, 'Value'); -[list_len,~] = size(entries); -[filename,filepath,~]=uigetfile(fullfile(pwd,'*.bld'),'Import Existing Population'); -fID = fopen(fullfile(filepath,filename),'r'); -list = textscan(fID,'%s\n'); -fclose(fID); -new_entries = list{1}; - -% Account for Empty List -if list_len == 0 || (list_len == 1 && strcmp(entries(1,:),'List of bld Data Files')) - for m=1:size(new_entries) - if m == 1 - entries = [new_entries{m}]; - else - entries = [entries; new_entries{m}]; - end - end - value = length(new_entries); -% Ensure that duplicate data files are not added -else - for m=1:size(new_entries) - on_list = false; - for n=1:size(entries) - if strcmp(entries(n,:),new_entries{m}) - on_list = true; - end - end - if ~on_list - entries = [entries; new_entries{m}]; - value = value+1; - end - end -end - -set(handles.list,'String',entries,'Value',value); - - - -% --- Executes on button press in save_pop. -function save_pop_Callback(hObject, eventdata, handles) -% hObject handle to save_pop (see GCBO) -% eventdata reserved - to be defined in a future version of MATLAB -% handles structure with handles and user data (see GUIDATA) -% -% Save the current list of files to be used again later (bld file) -value = get(handles.list, 'Value'); -entries = get(handles.list, 'String'); -if value < 2 - warndlg('Cannot save populations with fewer than 2 data files. Please add more data files to the list.', 'Populatin Too Small'); -else - [filename,filepath] = uiputfile('population.bld', 'Save Population As'); - fID = fopen(fullfile(filepath,filename),'w'); - fprintf(fID,'%s\n',entries{:}); - fclose(fID); -end diff --git a/gamma_bomb_2/code/pop_parent_fraction.m b/gamma_bomb_2/code/pop_parent_fraction.m deleted file mode 100644 index 78009d8..0000000 --- a/gamma_bomb_2/code/pop_parent_fraction.m +++ /dev/null @@ -1,52 +0,0 @@ -% Created by: Tom Morin -% Date: 5/28/15 - -function complete = pop_parent_fraction(files, directory) -% Generate a parent fraction curve from data from a population of several -% participants -% INPUT: files = bld files containing parent fraction data points - -% Read in data files, create matrices, and plot all matrices on one graph -complete = false; -cat_data = zeros(0,2); -[num_files,~] = size(files); -figure -for m=1:num_files - disp(files{m}); - A = dlmread(files{m}, '\t', 1, 0); - cat_data = cat(1, cat_data, A); - hold on; - %plot(A(:,1),A(:,2),'*r'); -end -title('Parent Fraction from Population Data','FontSize', 10,'FontWeight', 'bold'); -xlabel('Time (s)', 'FontSize', 7); -ylabel('Parent Fraction', 'FontSize', 7); -cat_data = sortrows(cat_data); - -pop_pf_mx = handle_repeated_data(cat_data); -plot(pop_pf_mx(:,1),pop_pf_mx(:,2),'*b'); - -% Write final Parent Fraction data file to csv_dir from gamma bomb gui -print_list({'sample-time[seconds]';'parent-fraction[1/1]'}, pop_pf_mx, fullfile(directory,'parentfraction.bld'),4); -complete = true; -end - - -function pop_pf_mx = handle_repeated_data(cat_data) -% Prepare the data points for curve fitting by ensuring that the data is a -% one-to-one function (one y-value for each x-value) -% If there are multiple y-values for one x-value, average them together - -pop_pf_mx(1,:) = cat_data(1,:); -num_repeats = 0; -for m=2:length(cat_data)-1 - if pop_pf_mx(end,1) ~= cat_data(m,1) - num_repeats = 0; - pop_pf_mx(end+1,:) = cat_data(m,:); - else - num_repeats = num_repeats + 1; - pop_pf_mx(end,2) = (pop_pf_mx(end,2).*num_repeats) + cat_data(m,2); - pop_pf_mx(end,2) = pop_pf_mx(end,2) ./ (num_repeats+1); - end -end -end \ No newline at end of file diff --git a/gamma_bomb_2/code/print_list.m b/gamma_bomb_2/code/print_list.m deleted file mode 100644 index 5f980d4..0000000 --- a/gamma_bomb_2/code/print_list.m +++ /dev/null @@ -1,55 +0,0 @@ -function print_list(hdr,M,ffile,sf) -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%% SCRIPT TO PRINT PMOD READABLE TEXT-FILES -%% Written by Daniel Chonde -%% Req: dynamic pet data (*.i files) -%% Output: Patlak i file (brain.i) -%% Written: 03/16/2012 Updated: -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% This script prints out data to a pmod readable tab delimited file. It is -% called by gammareader, the function used to automatically convert blood -% data to single files. -% hdr = the headings on the final output file -% M = a matrix of data -% ffile = the desired name of the output file -% sf = the desired significant figures for data in output file -% -% example: -% print_list({'sample-time[seconds]';plasma[kBq/cc]},data,'sampledata.bld') -% -% - -%% Set output number format -if nargin<4, sf=1;end -pformat=['%.' num2str(sf) 'f']; - -%% Load Row Data -if iscolumn(hdr),hdr=hdr';end -if size(M,2)~=size(hdr,2), error('data and data header are different lengths.');end - -%% Sort rows to make sure they are in ascending order -M=sortrows(M,1); - -%% Write Format Strings -hdr_str=''; -mat_str=''; -for m=1:size(hdr,2) -hdr_str=[hdr_str '%s']; -mat_str=[mat_str pformat]; -if m~=size(hdr,2), - hdr_str=[hdr_str '\t']; - mat_str=[mat_str '\t']; -else - hdr_str=[hdr_str '\n']; - mat_str=[mat_str '\n']; -end -end - -%% Print Data to file -fid=fopen(ffile,'w'); -fprintf(fid,hdr_str,hdr{:}); -fprintf(fid,mat_str,M'); -fclose(fid); - - -end \ No newline at end of file diff --git a/gamma_bomb_2/code/read_i_hdr.m b/gamma_bomb_2/code/read_i_hdr.m deleted file mode 100644 index c474ce2..0000000 --- a/gamma_bomb_2/code/read_i_hdr.m +++ /dev/null @@ -1,51 +0,0 @@ -function V=read_i_hdr(fullf) - - -if nargin<1 -[file,path,~]=uigetfile('*.i.hdr', 'Pick one of the DICOM files'); -fullf=fullfile(path,file); -end - -[path,filename,ext]=fileparts(fullf); -switch ext - case '.i' -file=strcat(filename,ext,'.hdr'); - case '.hdr' -file=strcat(filename,ext); -% case '' -% file=strcat(filename,'i.hdr'); - otherwise -error('File does not have .i or .i.hdr') -end -if isempty(path) - path=pwd; -end - -fid = fopen(fullfile(path,file), 'r'); %Open file for reading -m=1; -while 1 -tline=fgetl(fid); -if tline==-1;break;end -oldhdr{m,1} = tline; -clear tline -m=m+1; -end -fclose(fid); -clear m - -V=struct; -spoldhdr=cell(1); -for m=1:size(oldhdr,1) - if size(regexp(oldhdr{m},':=','split'),2)>size(spoldhdr,2); spoldhdr{:,end+1}='';end - spoldhdr(m,:)=regexp(oldhdr{m},':=','split'); - spoldhdr{m,1}=strtrim(regexprep(spoldhdr{m,1},'(\<[a-z])','${upper($1)}')); - spoldhdr{m,1}=regexprep(spoldhdr{m,1},' ',''); - spoldhdr{m,1}=regexprep(spoldhdr{m,1},'[\!,[,],/,(,),\%,:]',''); -end -clear m - -for m=1:size(spoldhdr,1) - V.(spoldhdr{m,1})=strtrim(spoldhdr{m,2}); -end - -end \ No newline at end of file diff --git a/gamma_bomb_2/final_report.pdf b/gamma_bomb_2/final_report.pdf deleted file mode 100644 index f54d9e1..0000000 Binary files a/gamma_bomb_2/final_report.pdf and /dev/null differ