data: [40x19x110 double]
time: [1x40 double]
energy: [15 25 35 45 55 65 75 85 95 105 115 125 135 145 155 165 175 185 195]
camstr = ['Cam #',num2str(hxr.cam{ic}.camnum)];
%
chordsind = hxr.cam{ic}.chordsind;% indices for chords of camera # ic
%
hxrdata = hxr.exp.data(:,:,chordsind);
%
if ichord > 0,% display for one chord
hxrdata = hxrdata(:,:,ichord);
chordstr = ['Chord #',num2str(hxr.cam{ic}.chordsnum(ichord))];
nchords = 1;
else
chordsnum = hxr.cam{ic}.chordsnum;
nchords = length(chordsnum);
end
%
if ibin > 0,% display for one energy
nbins = 1;
hxrdata = hxrdata(:,ibin,:);
if ibin < length(hxr.exp.energy),
energystr = [num2str(hxr.exp.energy(ibin)),' < k (keV) < ',num2str(hxr.exp.energy(ibin + 1))];
else
energystr = ['k (keV) > ',num2str(hxr.exp.energy(ibin))];
end
else
energyind = find(hxr.exp.energy(1:end-1) >= hxr.param.kmin & hxr.exp.energy(2:end) <= hxr.param.kmax);
%
energy = hxr.exp.energy;
energysep = [3*energy(1) - energy(2),energy(2:end) + energy(1:end-1),3*energy(end) - energy(end-1)]/2;
denergy = diff(energysep);
nbins = length(energy);
%
energystr = [num2str(hxr.param.kmin),' < k (keV) < ',num2str(hxr.param.kmax)];
end0001 function [nobj,graphids] = iluke_hxrbremdisp(shotimes,sdiags,ikphot,irho,itheta,idir,axs,select) 0002 % 0003 % data: [40x19x110 double] 0004 % time: [1x40 double] 0005 % energy: [15 25 35 45 55 65 75 85 95 105 115 125 135 145 155 165 175 185 195] 0006 % 0007 % camstr = ['Cam #',num2str(hxr.cam{ic}.camnum)]; 0008 % % 0009 % chordsind = hxr.cam{ic}.chordsind;% indices for chords of camera # ic 0010 % % 0011 % hxrdata = hxr.exp.data(:,:,chordsind); 0012 % % 0013 % if ichord > 0,% display for one chord 0014 % hxrdata = hxrdata(:,:,ichord); 0015 % chordstr = ['Chord #',num2str(hxr.cam{ic}.chordsnum(ichord))]; 0016 % nchords = 1; 0017 % else 0018 % chordsnum = hxr.cam{ic}.chordsnum; 0019 % nchords = length(chordsnum); 0020 % end 0021 % % 0022 % if ibin > 0,% display for one energy 0023 % nbins = 1; 0024 % hxrdata = hxrdata(:,ibin,:); 0025 % if ibin < length(hxr.exp.energy), 0026 % energystr = [num2str(hxr.exp.energy(ibin)),' < k (keV) < ',num2str(hxr.exp.energy(ibin + 1))]; 0027 % else 0028 % energystr = ['k (keV) > ',num2str(hxr.exp.energy(ibin))]; 0029 % end 0030 % else 0031 % energyind = find(hxr.exp.energy(1:end-1) >= hxr.param.kmin & hxr.exp.energy(2:end) <= hxr.param.kmax); 0032 % % 0033 % energy = hxr.exp.energy; 0034 % energysep = [3*energy(1) - energy(2),energy(2:end) + energy(1:end-1),3*energy(end) - energy(end-1)]/2; 0035 % denergy = diff(energysep); 0036 % nbins = length(energy); 0037 % % 0038 % energystr = [num2str(hxr.param.kmin),' < k (keV) < ',num2str(hxr.param.kmax)]; 0039 % end 0040 % 0041 nshotimes = length(shotimes); 0042 % 0043 nobj = [0,0,0]; 0044 graphids = {'','',''}; 0045 % 0046 if nshotimes > 1,% display time dependence 0047 % 0048 % timesind = []; 0049 % shotimesmask = false(1,nshotimes); 0050 % timesep = [3*hxr.exp.time(1) - hxr.exp.time(2),hxr.exp.time(2:end) + hxr.exp.time(1:end-1),3*hxr.exp.time(end) - hxr.exp.time(end-1)]/2; 0051 % dt = diff(timesep); 0052 % xlim = [timesep(1),timesep(end)]; 0053 % % 0054 % for ishotime = 1:nshotimes, 0055 % it = find(hxr.exp.time == shotimes(ishotime)); 0056 % % 0057 % if length(it) == 1, 0058 % timesind = [timesind,it]; 0059 % shotimesmask(ishotime) = true; 0060 % end 0061 % end 0062 % % 0063 % shotimes = shotimes(shotimesmask); 0064 % % 0065 % hxrdata = hxrdata(timesind,:,:); 0066 % dt = dt(timesind); 0067 % % 0068 % nt = length(timesind); 0069 % % 0070 % if ichord == 0,% display for all chords 0071 % % 0072 % % TODO 0073 % % 0074 % else% display for one chord 0075 % % 0076 % if ibin == 0,% display for all energies - show global count rate and photon temperature 0077 % % 0078 % Nphotons = repmat(dt,[nbins,1]).*hxrdata.';% number of photons measured per energy bin 0079 % % 0080 % hxrsum = sum(hxrdata(:,energyind),2).'; 0081 % ehxrsum = sqrt(hxrsum./dt); 0082 % % 0083 % [Tphexp,eTphexp] = tphfit_dke_yp(repmat(energy(energyind).',[1,nt]),Nphotons(energyind,:),sqrt(Nphotons(energyind,:))); 0084 % Tphmax = 2*Tphexp(find(hxrsum == max(hxrsum),1,'first')); 0085 % Tphmask = Tphexp > 0 & Tphexp < Tphmax & eTphexp < Tphexp; 0086 % Tphexp(~Tphmask) = NaN; 0087 % eTphexp(~Tphmask) = NaN; 0088 % % 0089 % graph1D_jd(shotimes,Tphexp,0,0,'t (s)','T_{ph.} (keV)','',NaN,xlim,'','-','o','r',2,style,axs(1),1,NaN,NaN,eTphexp);%Poloidal plot 0090 % ylim = get(axs(1),'ylim'); 0091 % if ~isnan(Tphmax), 0092 % set(axs(1),'ylim',[0,min([ylim(2),Tphmax])]);drawnow 0093 % end 0094 % % 0095 % nobj(1) = 1; 0096 % graphids{1} = ['Tph - ',camstr,'; ',chordstr,'; ',energystr]; 0097 % % 0098 % graph1D_jd(shotimes,hxrsum,0,0,'t (s)','CR (Counts/s)','',NaN,xlim,'','-','o','r',2,style,axs(2),1,NaN,NaN,ehxrsum);%Poloidal plot 0099 % % 0100 % nobj(2) = 1; 0101 % graphids{2} = ['CR - ',camstr,'; ',chordstr,'; ',energystr]; 0102 % % 0103 % else% display for one energy - display specific count rate 0104 % % 0105 % ehxrdata = sqrt(hxrdata./dt.'); 0106 % % 0107 % graph1D_jd(shotimes,hxrdata,0,0,'t (s)','CR (Counts/s)','',NaN,xlim,'','-','o','r',2,style,axs(1),1,NaN,NaN,ehxrdata);%Poloidal plot 0108 % % 0109 % nobj(1) = 1; 0110 % graphids{1} = ['CR - ',camstr,'; ',chordstr,'; ',energystr]; 0111 % % 0112 % end 0113 % % 0114 % end 0115 % 0116 else% display for one time step 0117 % 0118 timestr = ['t = ',num2str(shotimes,'%6.4f')]; 0119 Zbremplasma = sdiags{1}.Zbremplasma; 0120 equilHXR = sdiags{1}.equilHXR; 0121 % 0122 [nrho,ntheta] = size(Zbremplasma.theta_hxr); 0123 % 0124 bremdata = {Zbremplasma.brem_forward,Zbremplasma.brem_backward,Zbremplasma.brem_perp,Zbremplasma.brem_4pi}; 0125 % 0126 if ikphot == 0, 0127 % 0128 0129 % 0130 else 0131 % 0132 kphotstr = ['Photon energy k = ',num2str(Zbremplasma.kphot(ikphot)),' keV']; 0133 % 0134 if idir == 0, 0135 % 0136 0137 0138 % 0139 else% plot 2D FEB for a given energy and direction 0140 % 0141 dirstr = [select.dir{1 + idir},' direction']; 0142 % 0143 if irho == 0 && itheta == 0, 0144 % 0145 bremdisp = NaN(nrho,ntheta+1); 0146 % 0147 for ir = 1:nrho, 0148 for it = 1:ntheta, 0149 bremdisp(ir,it) = bremdata{idir}{ir,it}(ikphot); 0150 end 0151 end 0152 % 0153 bremdisp(:,end) = bremdisp(:,1);% for better plot 0154 % 0155 graph2D_jd(equilHXR.Xx,equilHXR.Xy,bremdisp.','x','y','',NaN,NaN,1,NaN,NaN,0,'-',NaN,0.5,select.style,1,NaN,axs(1)); 0156 axis(axs(1),'equal'); 0157 % 0158 nobj(1) = 1; 0159 graphids{1} = ['FEB - 2D ',timestr,'; ',dirstr,'; ',kphotstr]; 0160 % 0161 elseif irho == 0 && itheta > 0, 0162 % 0163 thetastr = ['theta = ',num2str(equilHXR.mtheta(itheta))*180/pi,' deg.']; 0164 % 0165 bremdisp = NaN(1,nrho); 0166 % 0167 for ir = 1:nrho, 0168 bremdisp(ir) = bremdata{idir}{ir,itheta}(ikphot); 0169 end 0170 % 0171 graph1D_jd(equilHXR.xrhoT,bremdisp,0,0,'r/a','emissivity (ph. s^{-1} keV^{-1} st^{-1})','',NaN,NaN,NaN,'-','o','r',2,select.style,axs(1)); 0172 % 0173 nobj(1) = 1; 0174 graphids{1} = ['FEB - rho ',timestr,'; ',dirstr,'; ',kphotstr,'; ',thetastr]; 0175 % 0176 elseif irho > 0 && itheta == 0, 0177 % 0178 rhostr = ['rho = ',num2str(equilHXR.xrhoT)]; 0179 % 0180 bremdisp = NaN(1,ntheta+1); 0181 % 0182 for it = 1:ntheta, 0183 bremdisp(it) = bremdata{idir}{irho,it}(ikphot); 0184 end 0185 % 0186 bremdisp(end) = bremdisp(1);% for better plot 0187 % 0188 graph1D_jd(equilHXR.mtheta,bremdisp,0,0,'\theta (deg.)','emissivity (ph. s^{-1} keV^{-1} st^{-1})','',NaN,NaN,NaN,'-','o','r',2,select.style,axs(1)); 0189 % 0190 nobj(1) = 1; 0191 graphids{1} = ['FEB - theta ',timestr,'; ',dirstr,'; ',kphotstr,'; ',rhostr]; 0192 % 0193 end 0194 % 0195 end 0196 % 0197 end 0198 0199 0200 % hxrdata = reshape(hxrdata(it,:,:),[nbins,nchords]); 0201 % % 0202 % CRmin = min(hxrdata(hxrdata(:) > 0)); 0203 % % 0204 % dt = 1/CRmin;% assuming a single photon count occured somewhere... to improve from database! 0205 % % 0206 % Nphotons = dt*hxrdata;% number of photons measured per energy bin 0207 % % 0208 % if ichord == 0,% display for all chords 0209 % % 0210 % xlim = [chordsnum(1) - 1,chordsnum(end)+1]; 0211 % % 0212 % if ibin == 0,% display for all energies - show global count rate and photon temperature 0213 % % 0214 % hxrsum = sum(hxrdata(energyind,:),1); 0215 % % 0216 % ehxrsum = sqrt(hxrsum/dt); 0217 % % 0218 % [Tphexp,eTphexp] = tphfit_dke_yp(repmat(energy(energyind).',[1,nchords]),Nphotons(energyind,:),sqrt(Nphotons(energyind,:))); 0219 % Tphmax = 2*Tphexp(find(hxrsum == max(hxrsum),1,'first')); 0220 % Tphmask = Tphexp > 0 & Tphexp < Tphmax & eTphexp < Tphexp; 0221 % Tphexp(~Tphmask) = NaN; 0222 % eTphexp(~Tphmask) = NaN; 0223 % % 0224 % graph1D_jd(chordsnum,Tphexp,0,0,'chord #','T_{ph.} (keV)','',NaN,xlim,'','-','o','r',2,style,axs(1),1,NaN,NaN,eTphexp);%Poloidal plot 0225 % ylim = get(axs(1),'ylim'); 0226 % if ~isnan(Tphmax), 0227 % set(axs(1),'ylim',[0,min([ylim(2),Tphmax])]);drawnow 0228 % end 0229 % % 0230 % nobj(1) = 1; 0231 % graphids{1} = ['Tph - ',camstr,'; ',timestr,'; ',energystr]; 0232 % % 0233 % graph1D_jd(chordsnum,hxrsum,0,0,'chord #','CR (Counts/s)','',NaN,xlim,'','-','o','r',2,style,axs(2),1,NaN,NaN,ehxrsum);%Poloidal plot 0234 % % 0235 % nobj(2) = 1; 0236 % graphids{2} = ['CR - ',camstr,'; ',timestr,'; ',energystr]; 0237 % % 0238 % else% display for one energy - display specific count rate 0239 % % 0240 % ehxrdata = sqrt(hxrdata/dt); 0241 % % 0242 % graph1D_jd(chordsnum,hxrdata,0,0,'chord #','CR (Counts/s)','',NaN,xlim,'','-','o','r',2,style,axs(1),1,NaN,NaN,ehxrdata);%Poloidal plot 0243 % % 0244 % nobj(1) = 1; 0245 % graphids{1} = ['CR - ',camstr,'; ',timestr,'; ',energystr]; 0246 % % 0247 % end 0248 % % 0249 % else% display spectrum 0250 % % 0251 % if ibin == 0,% display for all energies - show global count rate and photon temperature 0252 % % 0253 % xlim = [0,max(energysep)]; 0254 % % 0255 % [Tphexp,~,Aexp] = tphfit_dke_yp(energy(energyind).',Nphotons(energyind),sqrt(Nphotons(energyind))); 0256 % % 0257 % %log(Ei*Ni/dE) = a+b*Ei -> Tph = -1/b 0258 % % 0259 % ehxrdata = sqrt(hxrdata/dt); 0260 % % 0261 % hxrdatafit = Aexp*exp(-energy/Tphexp).*denergy./energy/dt; 0262 % % 0263 % graph1D_jd(energy,hxrdata,0,1,'k (keV)','CR (Counts/s)','',NaN,xlim,'','-','o','r',2,style,axs(1),1,NaN,NaN,ehxrdata);%Poloidal plot 0264 % graph1D_jd(energy,hxrdatafit,0,1,'','','',NaN,xlim,'','-','none','b',1,style,axs(1));%Poloidal plot 0265 % graph1D_jd(hxr.param.kmin*[1,1],get(axs(1),'ylim'),0,1,'','','',NaN,xlim,'','--','none','k',2,style,axs(1));%Poloidal plot 0266 % graph1D_jd(hxr.param.kmax*[1,1],get(axs(1),'ylim'),0,1,'','','',NaN,xlim,'','--','none','k',2,style,axs(1));%Poloidal plot 0267 % % 0268 % nobj(1) = 1; 0269 % graphids{1} = ['dCRdE - ',camstr,'; ',timestr,'; ',chordstr]; 0270 % % 0271 % else% display for one energy - display specific count rate 0272 % % 0273 % % TODO 0274 % % 0275 % end 0276 % % 0277 % end 0278 end 0279 %