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 mod_colors
21 : use colors_def
22 : use math_lib
23 : use const_def, only: dp, mesa_data_dir
24 : use utils_lib
25 :
26 : implicit none
27 : private
28 : public :: do_colors_init, free_colors_all, Eval_Colors
29 :
30 : contains
31 :
32 2 : subroutine do_colors_init(num_files,fnames,num_colors,ierr)
33 : integer,intent(in) :: num_files
34 : integer,dimension(:),intent(in) :: num_colors
35 : character(len=*),dimension(:),intent(in) :: fnames
36 : character(len=strlen) :: fname
37 : type (lgt_list), pointer :: thead =>null()
38 :
39 : integer, intent(out) :: ierr
40 : integer :: i
41 :
42 2 : ierr = 0
43 2 : num_thead=num_files
44 :
45 2 : if(num_thead<1)THEN
46 0 : color_is_initialized=.false.
47 : !Not a failure just dont have any files to read
48 : ierr=0
49 0 : return
50 : END IF
51 :
52 4 : ALLOCATE(thead_all(1:num_thead),STAT=ierr)
53 2 : if (ierr /= 0) then
54 0 : write(*,*) 'colors can not allocate memory in do_colors_init'
55 0 : ierr=-1
56 0 : return
57 : end if
58 :
59 2 : bc_total_num_colors=0
60 4 : do i=1,num_thead
61 2 : allocate(thead_all(i)%thead)
62 :
63 2 : thead=>thead_all(i)%thead
64 :
65 2 : fname=fnames(i)
66 2 : thead_all(i)%n_colors=num_colors(i)
67 : if(len(fname)==0)THEN
68 : exit
69 : end if
70 :
71 2 : if(thead_all(i)%n_colors<1)THEN
72 0 : write(*,*) "num_colors must be > 0"
73 0 : ierr=-1
74 0 : return
75 : end if
76 :
77 2 : call init_colors(fname,thead,thead_all(i)%color_names,thead_all(i)%n_colors,ierr)
78 :
79 2 : thead_all(i)%thead=>thead
80 2 : bc_total_num_colors=bc_total_num_colors+thead_all(i)%n_colors
81 4 : nullify(thead)
82 : end do
83 :
84 2 : color_is_initialized=.true.
85 :
86 : end subroutine do_colors_init
87 :
88 :
89 2 : subroutine init_colors(fname, thead, col_names, n_colors, ierr)
90 : integer, intent(out) :: ierr
91 : character(len=*),intent(in) :: fname
92 : type (lgt_list), pointer :: thead
93 : character(len=*),dimension(:) :: col_names
94 : integer, intent(in) :: n_colors
95 2 : call Read_Colors_Data(fname, thead, col_names, n_colors, ierr)
96 2 : end subroutine init_colors
97 :
98 2 : subroutine free_colors_all
99 : integer :: i
100 : type (lgt_list), pointer :: thead => null()
101 :
102 4 : do i=1,num_thead
103 4 : if(associated(thead_all(i)%thead))then
104 2 : thead=>thead_all(i)%thead
105 2 : call free_colors(thead)
106 : end if
107 : end do
108 2 : deallocate(thead_all)
109 :
110 2 : end subroutine free_colors_all
111 :
112 2 : subroutine free_colors(thead)
113 : type (lgt_list), pointer :: thead
114 : type(lgt_list),pointer :: tlist => null(), tnxt => null()
115 :
116 2 : tlist => thead
117 138 : do while (associated(tlist))
118 136 : tnxt => tlist% nxt
119 136 : call free_glist(tlist% glist)
120 136 : deallocate(tlist)
121 136 : tlist => tnxt
122 : end do
123 2 : nullify(thead)
124 :
125 : contains
126 :
127 136 : subroutine free_glist(gptr)
128 : type (lgg_list), pointer :: gptr
129 : type (lgg_list), pointer :: glist => null()
130 : type (lgg_list), pointer :: gnxt => null()
131 136 : glist => gptr
132 1076 : do while (associated(glist))
133 940 : gnxt => glist% nxt
134 940 : call free_zlist(glist% zlist)
135 940 : deallocate(glist)
136 940 : glist => gnxt
137 : end do
138 136 : end subroutine free_glist
139 :
140 940 : subroutine free_zlist(zptr)
141 : type (lgz_list), pointer :: zptr
142 : type (lgz_list), pointer :: zlist => null()
143 : type (lgz_list), pointer :: znxt => null()
144 940 : zlist => zptr
145 17512 : do while (associated(zlist))
146 16572 : znxt => zlist% nxt
147 16572 : deallocate(zlist)
148 16572 : zlist => znxt
149 : end do
150 940 : end subroutine free_zlist
151 :
152 : end subroutine free_colors
153 :
154 :
155 : subroutine show_tree(thead)
156 : type (lgt_list), pointer :: thead
157 : type (lgt_list), pointer :: tlist => null(), tnxt => null()
158 :
159 : tlist => thead
160 : do while (associated(tlist))
161 : write(*,*) 'tlist lgt', exp10(tlist% lgt)
162 : tnxt => tlist% nxt
163 : call show_glist(tlist% glist)
164 : tlist => tnxt
165 : end do
166 :
167 : contains
168 :
169 : subroutine show_glist(gptr)
170 : type (lgg_list), pointer :: gptr
171 : type (lgg_list), pointer :: glist => null(), gnxt => null()
172 :
173 : glist => gptr
174 : do while (associated(glist))
175 : write(*,*) 'glist% lgg', glist% lgg
176 : gnxt => glist% nxt
177 : call show_zlist(glist% zlist)
178 : glist => gnxt
179 : end do
180 : end subroutine show_glist
181 :
182 : subroutine show_zlist(zptr)
183 : type (lgz_list), pointer :: zptr
184 : type (lgz_list), pointer :: zlist => null()
185 : type (lgz_list), pointer :: znxt => null()
186 : zlist => zptr
187 : do while (associated(zlist))
188 : write(*,*) 'zlist% lgz', zlist% lgz
189 : znxt => zlist% nxt
190 : zlist => znxt
191 : end do
192 : end subroutine show_zlist
193 :
194 : end subroutine show_tree
195 :
196 :
197 2 : subroutine Read_One_Colors_Data(fname, thead, n_colors, col_names, ierr)
198 : integer, intent(out) :: ierr ! 0 means ok
199 : type (lgt_list), pointer,intent(inout) :: thead
200 :
201 : ! read file and build lists
202 : integer :: ios, cnt
203 : integer,intent(in) :: n_colors
204 : character (len=*) :: fname
205 2 : real(dp) :: lgt, lgg, lgz
206 : type (lgg_list), pointer :: glist => null()
207 : type (lgt_list), pointer :: tlist => null()
208 : type (lgz_list), pointer :: zlist => null()
209 42 : real(dp), dimension(max_num_bcs_per_file) :: colors
210 : character(len=*),dimension(:),intent(out) :: col_names
211 4 : character(len=256) :: tmp_cols(3+n_colors)
212 :
213 : character(len=4096) :: tmp
214 : integer :: num_entries, num_made, IO_UBV
215 :
216 : include 'formats'
217 :
218 : ! Try local folder first
219 2 : open(NEWUNIT=IO_UBV, FILE=trim(fname), ACTION='READ', STATUS='OLD', IOSTAT=ios)
220 2 : if(ios/=0) THEN
221 : ! Try colors data dir
222 1 : open(NEWUNIT=IO_UBV, FILE=trim(mesa_data_dir)//'/colors_data/'//trim(fname), ACTION='READ', STATUS='OLD', IOSTAT=ios)
223 1 : if (ios /= 0) then ! failed to open lcb98cor.dat
224 0 : write(*,*) 'colors_init: failed to open ' // trim(fname)
225 0 : write(*,*) 'or ',trim(mesa_data_dir)//'/colors_data/'//trim(fname)
226 0 : write(*,*)
227 0 : write(*,*) 'Please check that the file exists in your local work directory'
228 0 : write(*,*) 'or in ',trim(mesa_data_dir)//'/colors_data/'
229 0 : ierr = 1; return
230 : end if
231 : end if
232 :
233 2 : ierr = 0
234 2 : num_entries = 0
235 2 : cnt = 0
236 2 : tmp = ''
237 : !First line should be a header and containing the name of the colours
238 : !#teff logg m_div_h col_1 col_2 etc
239 2 : read(IO_UBV,'(a)') tmp
240 :
241 2 : call split_line(tmp,3+n_colors,tmp_cols)
242 :
243 24 : col_names(1:n_colors) = tmp_cols(4:n_colors+3)
244 :
245 49716 : do while (.true.)
246 16574 : read(IO_UBV,fmt=*,iostat=ios) lgt, lgg, lgz, colors(1:n_colors)
247 16574 : if (ios /= 0) exit
248 16572 : cnt = cnt + 1
249 16572 : lgt = log10(lgt)
250 :
251 16572 : if(cnt==1) thead%lgt = lgt
252 :
253 16572 : call get_tlist(thead, lgt, tlist, ierr)
254 16572 : if (ierr /= 0) exit
255 :
256 16572 : call get_glist(tlist% glist, lgg, glist, ierr)
257 16572 : if (ierr /= 0) exit
258 :
259 16572 : call get_zlist(glist% zlist, lgz, zlist, num_entries, ierr)
260 16572 : if (ierr /= 0) exit
261 :
262 16572 : if(zlist% colors(1) > -1d98) then
263 0 : write(*,*) "Warning found duplicated color data for (T, log g, M/H)=", 10**lgT, lgg, lgz
264 : end if
265 :
266 348014 : zlist% colors = colors
267 :
268 : end do
269 :
270 2 : close(IO_UBV)
271 2 : if (ierr /= 0) return
272 :
273 2 : num_made = 0
274 :
275 2 : tlist => thead
276 2 : lgt = 1d99
277 138 : do while (associated(tlist))
278 136 : if (tlist% lgt >= lgt) then ! bad tlist order
279 0 : ierr = -1; return
280 : end if
281 136 : lgt = tlist% lgt
282 136 : glist => tlist% glist
283 136 : lgg = 1d99
284 1076 : do while (associated(glist))
285 940 : if (glist% lgg >= lgg) then ! bad glist order
286 0 : ierr = -2; return
287 : end if
288 940 : lgg = glist% lgg
289 940 : zlist => glist% zlist
290 940 : lgz = 1d99
291 17512 : do while (associated(zlist))
292 16572 : if (zlist% lgz >= lgz) then ! bad zlist order
293 0 : ierr = -3; return
294 : end if
295 16572 : lgz = zlist% lgz
296 16572 : num_made = num_made + 1
297 16572 : zlist => zlist% nxt
298 : end do
299 940 : glist => glist% nxt
300 : end do
301 136 : tlist => tlist% nxt
302 : end do
303 :
304 2 : if(num_entries /= num_made)then
305 0 : write(0,*) "Error fond less colors than expected ",num_entries,num_made
306 0 : stop
307 : end if
308 :
309 :
310 : end subroutine Read_One_Colors_Data
311 :
312 :
313 2 : subroutine Read_Colors_Data(fname, thead, col_names, n_colors, ierr)
314 : integer, intent(out) :: ierr ! 0 means ok
315 : type (lgt_list), pointer,intent(inout) :: thead
316 : character (len=*),intent(in) :: fname
317 : character(len=*),dimension(:),intent(out) :: col_names
318 : integer, intent(in) :: n_colors
319 :
320 2 : Call Read_One_Colors_Data(fname, thead, n_colors, col_names, ierr)
321 2 : if (ierr /= 0) THEN
322 0 : write(*,*) "Read_Colors_One_Data error"
323 0 : stop
324 : end if
325 :
326 2 : end subroutine Read_Colors_Data
327 :
328 16572 : subroutine get_tlist(head, lgt, tlist, ierr)
329 : type (lgt_list), pointer :: head
330 : real(dp), intent(in) :: lgt
331 : type (lgt_list), pointer :: tlist
332 : integer, intent(out) :: ierr ! 0 means ok
333 :
334 : type (lgt_list), pointer :: t1=>null(), t2=>null()
335 :
336 16572 : ierr = 0
337 :
338 16572 : if (.not. associated(head)) then ! first time
339 0 : if (.not. alloc_tlist()) return
340 0 : head => tlist
341 0 : return
342 : end if
343 :
344 16572 : if (head% lgt == lgt) then ! matches head of list
345 804 : tlist => head
346 804 : return
347 : end if
348 :
349 15768 : if (head% lgt < lgt) then ! becomes new head of list
350 134 : if (.not. alloc_tlist()) return
351 134 : tlist% nxt => head
352 134 : head => tlist
353 134 : return
354 : end if
355 :
356 : ! check list
357 15634 : t1 => head
358 627696 : do while (associated(t1% nxt))
359 627696 : t2 => t1% nxt
360 627696 : if (t2% lgt == lgt) then
361 15634 : tlist => t2; return
362 : end if
363 612062 : if (t2% lgt < lgt) then ! insert new one before t2
364 0 : if (.not. alloc_tlist()) return
365 0 : tlist% nxt => t2
366 0 : t1% nxt => tlist
367 0 : return
368 : end if
369 612062 : t1 => t2
370 : end do
371 : ! add to end of list after t1
372 0 : if (.not. alloc_tlist()) return
373 0 : t1% nxt => tlist
374 :
375 : contains
376 :
377 134 : logical function alloc_tlist()
378 : integer :: istat
379 134 : allocate(tlist,stat=istat)
380 134 : if (istat /= 0) then
381 0 : alloc_tlist = .false.; ierr = -1; return
382 : end if
383 134 : nullify(tlist% glist)
384 134 : nullify(tlist% nxt)
385 134 : tlist% lgt = lgt
386 134 : alloc_tlist = .true.
387 134 : end function alloc_tlist
388 :
389 : end subroutine get_tlist
390 :
391 :
392 16572 : subroutine get_glist(head, lgg, glist, ierr)
393 : type (lgg_list), pointer :: head
394 : real(dp), intent(in) :: lgg
395 : type (lgg_list), pointer :: glist
396 : integer, intent(out) :: ierr
397 :
398 : type (lgg_list), pointer :: g1 => null(), g2 => null()
399 :
400 16572 : ierr = 0
401 :
402 16572 : if (.not. associated(head)) then ! first time
403 136 : if (.not. alloc_glist()) return
404 136 : head => glist
405 136 : return
406 : end if
407 :
408 16436 : if (head% lgg == lgg) then ! matches head of list
409 2316 : glist => head
410 2316 : return
411 : end if
412 :
413 14120 : if (head% lgg < lgg) then ! becomes new head of list
414 772 : if (.not. alloc_glist()) return
415 772 : glist% nxt => head
416 772 : head => glist
417 772 : return
418 : end if
419 :
420 : ! check list
421 13348 : g1 => head
422 55258 : do while (associated(g1% nxt))
423 55252 : g2 => g1% nxt
424 55252 : if (g2% lgg == lgg) then
425 13316 : glist => g2; return
426 : end if
427 41936 : if (g2% lgg < lgg) then ! insert new one before g2
428 26 : if (.not. alloc_glist()) return
429 26 : glist% nxt => g2
430 26 : g1% nxt => glist
431 26 : return
432 : end if
433 41910 : g1 => g2
434 : end do
435 : ! add to end of list after g1
436 6 : if (.not. alloc_glist()) return
437 6 : g1% nxt => glist
438 :
439 : contains
440 :
441 940 : logical function alloc_glist()
442 : integer :: istat
443 940 : allocate(glist,stat=istat)
444 940 : if (istat /= 0) then ! allocate failed in alloc_glist
445 0 : alloc_glist = .false.; ierr = -1; return
446 : end if
447 940 : nullify(glist% zlist)
448 940 : nullify(glist% nxt)
449 940 : glist% lgg = lgg
450 940 : alloc_glist = .true.
451 940 : end function alloc_glist
452 :
453 : end subroutine get_glist
454 :
455 16572 : subroutine get_zlist(head, lgz, zlist, num_entries, ierr)
456 : type (lgz_list), pointer :: head
457 : real(dp), intent(in) :: lgz
458 : type (lgz_list), pointer :: zlist
459 : integer, intent(out) :: ierr ! 0 means ok
460 : integer,intent(inout) :: num_entries
461 :
462 : type (lgz_list), pointer :: z1=>null(), z2=>null()
463 :
464 16572 : ierr = 0
465 :
466 16572 : if (.not. associated(head)) then ! first time
467 940 : if (.not. alloc_zlist()) return
468 940 : head => zlist
469 940 : return
470 : end if
471 :
472 15632 : if (head% lgz == lgz) then ! matches head of list
473 0 : zlist => head
474 0 : return
475 : end if
476 :
477 15632 : if (head% lgz < lgz) then ! becomes new head of list
478 5264 : if (.not. alloc_zlist()) return
479 5264 : zlist% nxt => head
480 5264 : head => zlist
481 5264 : return
482 : end if
483 :
484 : ! check list
485 10368 : z1 => head
486 65164 : do while (associated(z1% nxt))
487 54796 : z2 => z1% nxt
488 54796 : if (z2% lgz == lgz) then
489 0 : zlist => z2; return
490 : end if
491 54796 : if (z2% lgz < lgz) then ! insert new one before z2
492 0 : if (.not. alloc_zlist()) return
493 0 : zlist% nxt => z2
494 0 : z1% nxt => zlist
495 0 : return
496 : end if
497 54796 : z1 => z2
498 : end do
499 : ! add to end of list after z1
500 10368 : if (.not. alloc_zlist()) return
501 10368 : z1% nxt => zlist
502 :
503 : contains
504 :
505 16572 : logical function alloc_zlist()
506 : integer :: istat
507 348012 : allocate(zlist,stat=istat)
508 16572 : if (istat /= 0) then
509 0 : alloc_zlist = .false.; ierr = -1; return
510 : end if
511 16572 : nullify(zlist% nxt)
512 16572 : num_entries=num_entries+1
513 16572 : zlist% lgz = lgz
514 16572 : alloc_zlist = .true.
515 16572 : end function alloc_zlist
516 :
517 : end subroutine get_zlist
518 :
519 :
520 34 : subroutine Eval_Colors(log_Teff,log_g, M_div_h_in, results, thead, n_colors, ierr)
521 : real(dp), intent(in) :: log_Teff ! log10 of surface temp
522 : real(dp), intent(in) :: M_div_h_in ! [M/H]
523 : real(dp),dimension(:), intent(out) :: results
524 : real(dp), intent(in) :: log_g
525 : integer, intent(in) :: n_colors
526 : integer, intent(out) :: ierr
527 : type (lgt_list), pointer,intent(inout) :: thead
528 :
529 : !real(dp), parameter :: Zsol = 0.02d0, colors_bol_sun = 4.746d0
530 :
531 34 : real(dp) :: lgg, lgz, lgt, alfa, beta
532 1394 : real(dp),dimension(max_num_bcs_per_file) :: results1, results2
533 : type (lgt_list), pointer :: tlist => null(), tnxt => null()
534 :
535 34 : lgg=log_g
536 34 : lgz=M_div_h_in
537 34 : lgt=log_Teff
538 :
539 34 : if (.not. associated(thead)) then
540 33 : ierr = -1; return
541 : end if
542 :
543 34 : ierr = 0
544 :
545 : ! write(*,*) log_Teff,log_g, M_div_h_in
546 :
547 34 : if (lgt >= thead% lgt) then ! Error out
548 21 : results = -1d99
549 : return
550 : else
551 :
552 33 : tlist => thead
553 1700 : do while (associated(tlist% nxt))
554 1699 : tnxt => tlist% nxt
555 1699 : if (lgt == tnxt% lgt) then ! use tnxt
556 : call get_glist_results( &
557 0 : tnxt% glist, lgg, lgz, results1, n_colors, ierr)
558 0 : if (ierr /= 0) return
559 0 : results = results1
560 : return
561 : end if
562 1699 : if (lgt >= tnxt% lgt) then ! interpolate between tlist and tnxt
563 32 : call get_glist_results(tlist% glist, lgg, lgz, results1, n_colors, ierr)
564 32 : if (ierr /= 0) return
565 32 : call get_glist_results(tnxt% glist, lgg, lgz, results2, n_colors, ierr)
566 32 : if (ierr /= 0) return
567 702 : if(any(results1(1:n_colors)<-1d50) .or.any(results2(1:n_colors)<-1d50) ) then
568 63 : results = -1d99
569 : return
570 : end if
571 29 : alfa = (lgt - tnxt% lgt) / (tlist% lgt - tnxt% lgt)
572 29 : beta = 1.d0 - alfa
573 348 : results(1:n_colors) = alfa * results1(1:n_colors) + beta * results2(1:n_colors)
574 : return
575 : end if
576 1667 : tlist => tnxt
577 : end do
578 : end if
579 :
580 : ! Below all available log t's
581 21 : results = -1d99
582 :
583 :
584 : end subroutine Eval_Colors
585 :
586 :
587 64 : subroutine get_glist_results(gptr, lgg, lgz, results, n_colors, ierr)
588 : type (lgg_list), pointer :: gptr
589 : real(dp), intent(in) :: lgg, lgz
590 : real(dp),dimension(:), intent(out) :: results
591 : integer,intent(in) :: n_colors
592 : integer,intent(out) :: ierr
593 :
594 : type (lgg_list), pointer :: glist => null(), gnxt => null()
595 2624 : real(dp),dimension(max_num_bcs_per_file) :: results1, results2
596 64 : real(dp) :: alfa, beta
597 :
598 64 : glist => gptr
599 64 : ierr = 0
600 64 : if (.not. associated(glist)) then
601 0 : write(*,*) 'bad glist for get_glist_results'
602 0 : ierr = -1
603 60 : return
604 : end if
605 :
606 64 : if (lgg >= glist% lgg) then ! use the largest lgg
607 42 : results = -1d99
608 : return
609 : end if
610 :
611 160 : do while (associated(glist% nxt))
612 156 : gnxt => glist% nxt
613 156 : if (lgg == gnxt% lgg) then ! use gnxt
614 0 : call get_zlist_results(gnxt% zlist, lgz, results1, n_colors, ierr)
615 0 : if (ierr /= 0) return
616 0 : results=results1
617 : return
618 : end if
619 156 : if (lgg >= gnxt% lgg) then ! interpolate between lgg's
620 58 : call get_zlist_results(glist% zlist, lgz, results1, n_colors, ierr)
621 58 : if (ierr /= 0) return
622 58 : call get_zlist_results(gnxt% zlist, lgz, results2, n_colors, ierr)
623 58 : if (ierr /= 0) return
624 :
625 1392 : if(any(results1(1:n_colors)<-1d50) .or.any(results2(1:n_colors)<-1d50) ) then
626 0 : results = -1d99
627 : return
628 : end if
629 58 : alfa = (lgg - gnxt% lgg) / (glist% lgg - gnxt% lgg)
630 58 : beta = 1.d0 - alfa
631 696 : results(1:n_colors) = alfa * results1(1:n_colors) + beta * results2(1:n_colors)
632 : return
633 : end if
634 98 : glist => gnxt
635 : end do
636 :
637 : ! below all available log g's
638 84 : results = -1d99
639 :
640 : end subroutine get_glist_results
641 :
642 :
643 116 : subroutine get_zlist_results(zptr, lgz, results, n_colors, ierr)
644 : type (lgz_list), pointer :: zptr
645 : real(dp), intent(in) :: lgz
646 : real(dp),dimension(:), intent(out) :: results
647 : integer, intent(in) :: n_colors
648 : integer, intent(out) :: ierr
649 :
650 : type (lgz_list), pointer :: zlist => null(), znxt => null()
651 116 : real(dp) :: alfa, beta
652 :
653 116 : zlist => zptr
654 :
655 116 : ierr = 0
656 116 : if (.not. associated(zlist)) then
657 0 : write(*,*) 'bad zlist for get_zlist_results'
658 0 : ierr = -1
659 0 : return
660 : end if
661 :
662 116 : if (lgz >= zlist% lgz) then ! use the largest lgz
663 0 : results = -1d99
664 : return
665 : end if
666 :
667 580 : do while (associated(zlist% nxt))
668 580 : znxt => zlist% nxt
669 580 : if (lgz == znxt% lgz) then ! use znxt
670 1392 : results(1:n_colors) = znxt% colors(1:n_colors)
671 : return
672 : end if
673 464 : if (lgz >= znxt% lgz) then ! interpolate between zlist and znxt
674 0 : alfa = (lgz - znxt% lgz) / (zlist% lgz - znxt% lgz)
675 0 : beta = 1 - alfa
676 0 : results(1:n_colors) = alfa * zlist% colors(1:n_colors) + beta * znxt% colors(1:n_colors)
677 : return
678 : end if
679 464 : zlist => znxt
680 : end do
681 :
682 : ! below all available log z's
683 0 : results = -1d99
684 :
685 : end subroutine get_zlist_results
686 :
687 :
688 : end module mod_colors
689 :
|
