Line data Source code
1 : ! ***********************************************************************
2 : !
3 : ! Copyright (C) 2010-2019 The MESA Team
4 : !
5 : ! This program is free software: you can redistribute it and/or modify
6 : ! it under the terms of the GNU Lesser General Public License
7 : ! as published by the Free Software Foundation,
8 : ! either version 3 of the License, or (at your option) any later version.
9 : !
10 : ! This program is distributed in the hope that it will be useful,
11 : ! but WITHOUT ANY WARRANTY; without even the implied warranty of
12 : ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
13 : ! See the GNU Lesser General Public License for more details.
14 : !
15 : ! You should have received a copy of the GNU Lesser General Public License
16 : ! along with this program. If not, see <https://www.gnu.org/licenses/>.
17 : !
18 : ! ***********************************************************************
19 :
20 : module pgstar_power
21 :
22 : use star_private_def
23 : use const_def, only: dp
24 : use pgstar_support
25 : use star_pgstar
26 :
27 : implicit none
28 :
29 : contains
30 :
31 0 : subroutine power_plot(id, device_id, ierr)
32 : integer, intent(in) :: id, device_id
33 : integer, intent(out) :: ierr
34 :
35 : type (star_info), pointer :: s
36 :
37 : ierr = 0
38 0 : call get_star_ptr(id, s, ierr)
39 0 : if (ierr /= 0) return
40 :
41 0 : call pgslct(device_id)
42 0 : call pgbbuf()
43 0 : call pgeras()
44 :
45 : call do_power_plot(s, id, device_id, &
46 : s% pg% Power_xleft, s% pg% Power_xright, &
47 : s% pg% Power_ybot, s% pg% Power_ytop, .false., &
48 0 : s% pg% Power_title, s% pg% Power_txt_scale, ierr)
49 :
50 0 : call pgebuf()
51 :
52 : end subroutine power_plot
53 :
54 :
55 0 : subroutine do_power_plot(s, id, device_id, &
56 : winxmin, winxmax, winymin, winymax, subplot, &
57 : title, txt_scale, ierr)
58 : type (star_info), pointer :: s
59 : integer, intent(in) :: id, device_id
60 : real, intent(in) :: winxmin, winxmax, winymin, winymax
61 : logical, intent(in) :: subplot
62 : character (len=*), intent(in) :: title
63 : real, intent(in) :: txt_scale
64 : integer, intent(out) :: ierr
65 : call do_power_panel(s, id, device_id, &
66 : winxmin, winxmax, winymin, winymax, subplot, &
67 : title, txt_scale, s% pg% Power_xaxis_name, &
68 : s% pg% Power_xmin, s% pg% Power_xmax, s% pg% Power_xaxis_reversed, &
69 0 : s% pg% Power_ymin, s% pg% Power_ymax, .false., .true., ierr)
70 0 : end subroutine do_power_plot
71 :
72 :
73 0 : subroutine do_power_panel(s, id, device_id, &
74 : winxmin, winxmax, winymin, winymax, subplot, title, txt_scale, &
75 : xaxis_name, xaxis_min, xaxis_max, xaxis_reversed, ymin_in, ymax_in, &
76 : panel_flag, xaxis_numeric_labels_flag, ierr)
77 : use utils_lib
78 : use chem_def
79 : use net_def
80 : use const_def, only: Msun, Rsun
81 : use pgstar_colors
82 :
83 : type (star_info), pointer :: s
84 : integer, intent(in) :: id, device_id
85 : real, intent(in) :: &
86 : winxmin, winxmax, winymin, winymax, xaxis_min, xaxis_max, &
87 : ymin_in, ymax_in
88 : logical, intent(in) :: subplot
89 : character (len=*), intent(in) :: title, xaxis_name
90 : real, intent(in) :: txt_scale
91 : logical, intent(in) :: &
92 : xaxis_reversed, panel_flag, xaxis_numeric_labels_flag
93 : integer, intent(out) :: ierr
94 :
95 : real, allocatable, dimension(:) :: xvec, yvec
96 0 : real :: xmin, xmax, xleft, xright, dx, chScale, windy, &
97 0 : ymin, ymax, exp10_ymin, xmargin, &
98 : legend_xmin, legend_xmax, legend_ymin, legend_ymax
99 : integer :: lw, lw_sav, grid_min, grid_max, npts, nz
100 : integer, parameter :: num_colors = 20
101 : integer :: colors(num_colors)
102 :
103 : include 'formats'
104 0 : ierr = 0
105 0 : nz = s% nz
106 :
107 : colors(:) = [ &
108 : clr_MediumSlateBlue, clr_LightSkyBlue, clr_Goldenrod, clr_Lilac, &
109 : clr_Coral, clr_Crimson, clr_LightSkyGreen, clr_DarkGray, &
110 : clr_Tan, clr_IndianRed, clr_Gold, &
111 : clr_Teal, clr_Silver, clr_BrightBlue, clr_FireBrick, &
112 : clr_RoyalPurple, clr_SlateGray, clr_LightSteelBlue, &
113 0 : clr_Gray, clr_RoyalBlue ]
114 :
115 0 : chScale = txt_scale
116 :
117 0 : windy = winymax - winymin
118 :
119 0 : legend_xmin = winxmax
120 0 : legend_xmax = min(1.0, winxmax + (winxmax - winxmin)/5)
121 0 : legend_ymin = winymin
122 0 : legend_ymax = winymax
123 :
124 0 : allocate (xvec(nz), yvec(nz))
125 :
126 0 : xmargin = 0
127 : call set_xaxis_bounds(s, xaxis_name, xaxis_min, xaxis_max, &
128 : xaxis_reversed, xmargin, xvec, xmin, xmax, xleft, xright, dx, &
129 0 : grid_min, grid_max, npts, ierr)
130 :
131 0 : call pgsave
132 0 : call pgsch(txt_scale)
133 0 : if (ierr == 0) call plot(ierr)
134 0 : call pgunsa
135 :
136 0 : deallocate(xvec, yvec)
137 :
138 : contains
139 :
140 :
141 0 : subroutine plot(ierr)
142 0 : use rates_def
143 : integer, intent(out) :: ierr
144 :
145 : integer :: i, j, cnt
146 : logical, parameter :: dbg = .false.
147 0 : real(dp) :: max_power(num_categories), max_power_copy(num_categories)
148 0 : real :: ybot
149 :
150 : include 'formats'
151 :
152 0 : ymax = -1e9
153 0 : do i = 1, num_categories
154 0 : max_power(i) = maxval(s% eps_nuc_categories(i,grid_min:grid_max))
155 0 : if (max_power(i) > ymax) ymax = max_power(i)
156 0 : max_power_copy(i) = max_power(i)
157 : end do
158 :
159 0 : if (ymax < 1e-29) ymax = 1e-29
160 0 : ymax = log10(dble(ymax))
161 0 : if (ymax <= 4) then
162 0 : ymax = 4.3
163 0 : ymin = -4.1
164 : else
165 0 : ymax = ymax*1.1
166 0 : ymin = -0.1
167 : end if
168 :
169 0 : if (ymax_in /= -101) ymax = ymax_in
170 0 : if (ymin_in /= -101) ymin = ymin_in
171 :
172 0 : exp10_ymin = exp10(dble(ymin))
173 :
174 0 : lw = s% pg% pgstar_lw
175 0 : call pgqlw(lw_sav)
176 :
177 0 : call pgsave
178 0 : call pgsvp(legend_xmin, legend_xmax, legend_ymin, legend_ymax)
179 0 : call pgswin(0.0, 1.0, ymin, ymax)
180 0 : cnt = 0
181 0 : do j=1,num_categories
182 0 : i = maxloc(max_power(:),dim=1)
183 0 : cnt = power_line_legend(cnt,i)
184 0 : max_power(i) = -1d99
185 : end do
186 0 : call pgunsa
187 :
188 0 : call pgsvp(winxmin, winxmax, winymin, winymax)
189 0 : if (.not. panel_flag) then
190 0 : if (.not. subplot) then
191 0 : call show_model_number_pgstar(s)
192 0 : call show_age_pgstar(s)
193 : end if
194 0 : call show_title_pgstar(s, title)
195 0 : call pgsci(clr_Foreground)
196 0 : call show_xaxis_name(s,xaxis_name,ierr)
197 0 : if (ierr /= 0) return
198 : end if
199 :
200 0 : ybot = -0.05
201 0 : call pgswin(xleft, xright, ymin+ybot, ymax)
202 0 : call pgscf(1)
203 0 : call pgsci(clr_Foreground)
204 0 : if (xaxis_numeric_labels_flag) then
205 0 : call show_box_pgstar(s,'BCNST','BCNSTV')
206 : else
207 0 : call show_box_pgstar(s,'BCST','BCNSTV')
208 : end if
209 0 : call show_left_yaxis_label_pgstar(s,'log ergs/g/s')
210 :
211 0 : call pgslw(lw)
212 0 : cnt = 0
213 0 : do j=1,num_categories
214 0 : i = maxloc(max_power_copy(:),dim=1)
215 0 : cnt = power_line(cnt,i)
216 0 : max_power_copy(i) = -1d99
217 : end do
218 0 : call pgslw(lw_sav)
219 :
220 0 : if (.not. panel_flag) then ! show mix regions at bottom of plot
221 0 : call pgslw(10)
222 : call show_mix_regions_on_xaxis( &
223 0 : s,ymin+ybot,ymax,grid_min,grid_max,xvec)
224 : end if
225 :
226 0 : call show_pgstar_decorator(s%id,s% pg% power_use_decorator,s% pg% power_pgstar_decorator, 0, ierr)
227 :
228 :
229 0 : end subroutine plot
230 :
231 :
232 0 : integer function power_line(cnt, icat)
233 : integer, intent(in) :: cnt, icat
234 0 : real :: ymx
235 : integer :: iclr, k
236 0 : power_line = cnt
237 0 : ymx = maxval(s% eps_nuc_categories(icat,grid_min:grid_max))
238 0 : if (ymx < exp10_ymin) return
239 0 : iclr = modulo(icat,num_colors) +1
240 0 : power_line = cnt + 1
241 0 : call pgsci(colors(iclr))
242 0 : do k=1,nz
243 0 : yvec(k) = safe_log10(s% eps_nuc_categories(icat,k))
244 : end do
245 0 : call pgline(npts, xvec(grid_min:grid_max), yvec(grid_min:grid_max))
246 0 : end function power_line
247 :
248 :
249 0 : integer function power_line_legend(cnt, icat)
250 : integer, intent(in) :: cnt, icat
251 0 : real :: ymx, dx, dyline, ypos, xpts(2), ypts(2)
252 : integer :: iclr, num_max
253 0 : num_max = min(num_categories, s% pg% Power_legend_max_cnt)
254 0 : power_line_legend = cnt
255 0 : if (cnt >= num_max) return
256 0 : ymx = maxval(s% eps_nuc_categories(icat,grid_min:grid_max))
257 0 : if (ymx < exp10_ymin) return
258 0 : iclr = modulo(icat,num_colors) +1
259 0 : power_line_legend = cnt + 1
260 0 : call pgsci(colors(iclr))
261 0 : dx = 0.1
262 0 : dyline = (ymax-ymin)/num_max
263 0 : ypos = ymax - (cnt+0.5)*dyline
264 0 : xpts(1) = 1.3*dx
265 0 : xpts(2) = xpts(1) + 2.3*dx
266 0 : ypts = ypos + dyline*0.1
267 0 : call pgslw(lw)
268 0 : call pgline(2, xpts, ypts)
269 0 : call pgslw(lw_sav)
270 0 : call pgsci(clr_Foreground)
271 0 : call pgsch(txt_scale*s% pg% Power_legend_txt_scale_factor)
272 : call pgptxt(xpts(2) + dx, ypos, 0.0, 0.0, &
273 0 : trim(adjustl(category_name(icat))))
274 0 : power_line_legend = cnt + 1
275 0 : end function power_line_legend
276 :
277 :
278 : end subroutine do_power_panel
279 :
280 : end module pgstar_power
|
