3 f1 = figure(
'Units',
'inches',
'Position', [0.5 0.5 8.5 4], ...
4 'DefaultAxesFontName',
'times new roman',
'DefaultAxesFontWeight',
'bold',
'DefaultAxesFontSize', 7,
'DefaultTextFontWeight',
'bold',
'DefaultTextFontSize', 7,...
5 'PaperOrientation',
'portrait',
'PaperPositionMode',
'auto',
'DefaultLineLineWidth', 2);%,
'PaperPosition', [1 1 3.8 9]);%
14 subplot1(2, 3, 'Max', [0.97 0.95], 'Min', [0.05 0.1], 'Gap', [0.02 0.02]);
19 {
'../Execute/3D_Without_0.1hr v1.28/Output/',
' No mixed terms'},
20 {
'../Execute/3D_With_0.1hr v1.28/Output/',
' With mixed terms '},
24 {'../Execute/3D_Without_0.1hr v1.28/Output/', ' No mixed terms'},
25 {'../Execute/3D_With_0.1hr v1.28/Output/', ' With mixed terms '},
29 path1 = path_arr1{2,1};
30 path2 = path_arr1{1,1};
32 name1 = path_arr1{2,2};
33 name2 = path_arr1{1,2};
47 eval([
'load ''', [path1,
'Flux_E=', num2str(target_epc),
'_alpha', num2str(target_alpha)],
''''])
52 eval(['load
''', [path2,
'Flux_E=', num2str(target_epc),
'_alpha', num2str(target_alpha)],
''''])
57 if (max(size(Flux1) ~= size(Flux2)))
58 Flux2 = interp2(time2, L2, Flux2', time1
', L1, 'spline')';
61 [ND, AND, MND] = NormDiff(Flux1, Flux2, 2);
65 Flux_min = 10^c_min;%Flux_max - 10^delta_color;
66 Flux1(find(Flux1<Flux_min))=Flux_min;
67 log_Flux = log10(Flux1);
69 contourf(time_arr, L1, log_Flux
', n_colors);
75 text(max(xlim)-0.05, (min(ylim) + max(ylim))/2, 'Flux, log_{10}(#/sr/s/cm^2/KeV)
', 'Rotation
', 90, 'VerticalAlignment
', 'top
', 'HorizontalAlignment
', 'center
', 'FontWeight
', 'bold
');
76 %text(max(xlim), (min(ylim) + max(ylim))/2, '-3 \cdot 10^{-8} (c/MeV/cm)^3', 'Rotation', 90, 'VerticalAlignment', 'top', 'HorizontalAlignment', 'center', 'FontWeight', 'bold');
77 %xlabel(
'Time, day of Apr. 2001');
80 text(min(xlim),
max(ylim), [
' (D) '],
'VerticalAlignment',
'top',
'HorizontalAlignment',
'left',
'FontWeight',
'bold',
'Color',
'white');
81 if (target_alpha >= 88)
82 text(
max(xlim), min(ylim), [{name1}, {[' Energy = ', num2str(target_epc), ' MeV, equatorial pitch angle ']}],
'VerticalAlignment',
'bottom',
'HorizontalAlignment',
'right',
'FontWeight',
'bold',
'color', TextColor);
84 text(
max(xlim), min(ylim), [{name1}, {[' Energy = ', num2str(target_epc), ' MeV, Pitch angle = ', num2str(target_alpha), '^o ']}],
'VerticalAlignment',
'bottom',
'HorizontalAlignment',
'right',
'FontWeight',
'bold',
'color', TextColor);
91 Flux_min = 10^c_min;%Flux_max - 10^delta_color;
92 Flux2(find(Flux2<Flux_min))=Flux_min;
93 log_Flux = log10(Flux2);
95 contourf(time_arr, L1, log_Flux
', n_colors);
101 text(max(xlim)-0.05, (min(ylim) + max(ylim))/2, 'Flux, log_{10}(#/sr/s/cm^2/KeV)
', 'Rotation
', 90, 'VerticalAlignment
', 'top
', 'HorizontalAlignment
', 'center
', 'FontWeight
', 'bold
');
102 %text(max(xlim), (min(ylim) + max(ylim))/2, '-3 \cdot 10^{-8} (c/MeV/cm)^3', 'Rotation', 90, 'VerticalAlignment', 'top', 'HorizontalAlignment', 'center', 'FontWeight', 'bold');
103 %xlabel(
'Time, day of Apr. 2001');
104 xlabel(
'Time, days');
106 text(min(xlim),
max(ylim), [
' (E) '],
'VerticalAlignment',
'top',
'HorizontalAlignment',
'left',
'FontWeight',
'bold',
'Color',
'white');
107 if (target_alpha >= 88)
108 text(
max(xlim), min(ylim), [{name2}, {[' Energy = ', num2str(target_epc), ' MeV, equatorial pitch angle ']}],
'VerticalAlignment',
'bottom',
'HorizontalAlignment',
'right',
'FontWeight',
'bold',
'color', TextColor);
110 text(
max(xlim), min(ylim), [{name2}, {[' Energy = ', num2str(target_epc), ' MeV, Pitch angle = ', num2str(target_alpha), '^o ']}],
'VerticalAlignment',
'bottom',
'HorizontalAlignment',
'right',
'FontWeight',
'bold',
'color', TextColor);
116 contourf(time_arr, L1, ND
', n_colors);
122 text(max(xlim), (min(ylim) + max(ylim))/2, 'Normalized difference, %
', 'Rotation
', 90, 'VerticalAlignment
', 'top
', 'HorizontalAlignment
', 'center
', 'FontWeight
', 'bold
');
123 %text(max(xlim), (min(ylim) + max(ylim))/2, '-3 \cdot 10^{-8} (c/MeV/cm)^3', 'Rotation', 90, 'VerticalAlignment', 'top', 'HorizontalAlignment', 'center', 'FontWeight', 'bold');
124 %xlabel(
'Time, day of Apr. 2001');
125 xlabel(
'Time, days');
127 text(
max(xlim), min(ylim), {['(F_2 - F_1)/max_{at each t}(F_2/2 + F_1/2)']},
'VerticalAlignment',
'bottom',
'HorizontalAlignment',
'right',
'FontWeight',
'bold',
'color',
'black');
128 text(min(xlim),
max(ylim), [
' (F) '],
'VerticalAlignment',
'top',
'HorizontalAlignment',
'left',
'FontWeight',
'bold',
'Color',
'black');
130 %%save_picture(gca, [path1,
'../Plots/',
'ND_Flux_E=', num2str(target_epc),
'_alpha', num2str(target_alpha)],
'jpg',
'pdf',
'epsc',
'fig');
134 path1 = path_arr2{2,1};
135 path2 = path_arr2{1,1};
137 name1 = path_arr2{2,2};
138 name2 = path_arr2{1,2};
148 eval([
'load ''', path1,
'PSD_mu=', num2str(target_mu),
'_J', num2str(target_J),
''' -mat
'])
153 eval(['load
''', path2,
'PSD_mu=', num2str(target_mu),
'_J', num2str(target_J),
''' -mat
'])
158 if (max(size(PSD1) ~= size(PSD2)))
159 PSd2 = interp2(time2, L2, PSD2', time1
', L1, 'spline')';
162 [ND, AND, MND] = NormDiff(PSD1, PSD2, 2);
167 PSD1(find(PSD1<PSD_min))=PSD_min;
168 log_PSD = log10(PSD1);
170 contourf(time_arr, L1, log_PSD
', n_colors);
174 caxis([c_min c_max]);
176 text(max(xlim)-0.05, (min(ylim) + max(ylim))/2, 'PSD, log_{10}(3 \cdot 10^{-8} (c/cm/MeV)^3)', 'Rotation', 90, 'VerticalAlignment', 'top', 'HorizontalAlignment', 'center', 'FontWeight', 'bold');
177 %text(
max(xlim), (min(ylim) +
max(ylim))/2,
'-3 \cdot 10^{-8} (c/MeV/cm)^3',
'Rotation', 90,
'VerticalAlignment',
'top',
'HorizontalAlignment',
'center',
'FontWeight',
'bold');
178 %xlabel(
'Time, day of Apr. 2001');
179 xlabel(
'Time, days');
181 text(
max(xlim), min(ylim), [{name1}, {[' \mu = ', num2str(target_mu), ' MeV / G, J = ', num2str(target_J), 'MeV/c Re ']}],
'VerticalAlignment',
'bottom',
'HorizontalAlignment',
'right',
'FontWeight',
'bold',
'color', TextColor);
182 text(min(xlim),
max(ylim), [
' (A) '],
'VerticalAlignment',
'top',
'HorizontalAlignment',
'left',
'FontWeight',
'bold',
'Color',
'black');
186 PSD_min = 10^c_min;%PSD_max - 10^delta_color;
187 PSD2(find(PSD2<PSD_min))=PSD_min;
188 log_PSD = log10(PSD2);
190 contourf(time_arr, L1, log_PSD
', n_colors);
194 caxis([c_min c_max]);
196 text(max(xlim)-0.05, (min(ylim) + max(ylim))/2, 'PSD, log_{10}(3 \cdot 10^{-8} (c/cm/MeV)^3)', 'Rotation', 90, 'VerticalAlignment', 'top', 'HorizontalAlignment', 'center', 'FontWeight', 'bold');
197 %text(
max(xlim), (min(ylim) +
max(ylim))/2,
'-3 \cdot 10^{-8} (c/MeV/cm)^3',
'Rotation', 90,
'VerticalAlignment',
'top',
'HorizontalAlignment',
'center',
'FontWeight',
'bold');
198 %xlabel(
'Time, day of Apr. 2001');
199 xlabel(
'Time, days');
201 text(
max(xlim), min(ylim), [{name2}, {[' \mu = ', num2str(target_mu), ' MeV / G, J = ', num2str(target_J), 'MeV/c Re ']}],
'VerticalAlignment',
'bottom',
'HorizontalAlignment',
'right',
'FontWeight',
'bold',
'color', TextColor);
202 text(min(xlim),
max(ylim), [
' (B) '],
'VerticalAlignment',
'top',
'HorizontalAlignment',
'left',
'FontWeight',
'bold',
'Color',
'black');
206 contourf(time_arr, L1, ND
', n_colors);
212 text(max(xlim), (min(ylim) + max(ylim))/2, 'Normalized difference, %
', 'Rotation
', 90, 'VerticalAlignment
', 'top
', 'HorizontalAlignment
', 'center
', 'FontWeight
', 'bold
');
213 %text(max(xlim), (min(ylim) + max(ylim))/2, '-3 \cdot 10^{-8} (c/MeV/cm)^3', 'Rotation', 90, 'VerticalAlignment', 'top', 'HorizontalAlignment', 'center', 'FontWeight', 'bold');
214 %xlabel(
'Time, day of Apr. 2001');
215 xlabel(
'Time, days');
217 text(
max(xlim), min(ylim), {['(F_2 - F_1)/max_{at each t}(F_2/2 + F_1/2)']},
'VerticalAlignment',
'bottom',
'HorizontalAlignment',
'right',
'FontWeight',
'bold',
'color',
'black');
218 text(min(xlim),
max(ylim), [
' (C) '],
'VerticalAlignment',
'top',
'HorizontalAlignment',
'left',
'FontWeight',
'bold',
'Color',
'black');
220 %%save_picture(gca, [path1,
'../Plots/',
'ND_PSD_mu=', num2str(target_mu),
'J=', num2str(target_J)],
'jpg',
'pdf',
'epsc',
'fig');
222 save_picture(
gcf, [
'Figure_5'],
'jpg',
'pdf',
'epsc');
223 %save_picture(
gcf, [
'Figure_5_bw'],
'jpg',
'pdf',
'eps');
double max(double v1, double v2)
Return maximum.
void gcf(double *gammcf, double a, double x, double *gln)
Returns the incomplete gamma function Q(a, x) evaluated by its continued fraction representation as g...
Phase Space Density class.
void spline(double *x, double *y, int n, double yp1, double ypn, double *y2)