source: trunk/LMDZ.GENERIC/libf/phystd/ave_stelspec.F90 @ 603

Last change on this file since 603 was 374, checked in by emillour, 13 years ago

Generic GCM: Upgrade: The location of the 'datagcm' directory can now be given in the callphys.def file ( datadir = /absolute/path/to/datagcm ). Changed "datafile.h" into a F90 module "datafile_mod.F90" and spread this change to all routines that used to use "datafile.h".
EM

File size: 4.0 KB
Line 
1      subroutine ave_stelspec(STELLAR)
2
3!==================================================================
4!     
5!     Purpose
6!     -------
7!     Average the chosen high resolution stellar spectrum over the
8!     visible bands in the model.
9!     
10!     Authors
11!     -------
12!     Robin Wordsworth (2010).
13!
14!     Called by
15!     ---------
16!     setspv.F
17!     
18!     Calls
19!     -----
20!     none
21!     
22!==================================================================
23
24      use radinc_h, only: L_NSPECTV
25      use radcommon_h, only: BWNV, DWNV, tstellar
26      use datafile_mod, only: datadir
27
28      implicit none
29
30#include "callkeys.h"
31
32      real*8 STELLAR(L_NSPECTV)
33!      integer startype
34
35      integer Nfine
36      parameter(Nfine=5000)
37      integer ifine
38
39      real lam(Nfine)
40      real stel_f(Nfine)
41      real band
42      real dl
43
44      character(len=50)  :: file_id
45      character(len=100) :: file_path
46
47      real lam_temp
48      double precision stel_temp
49     
50      integer :: ios ! file opening/reading status
51
52      STELLAR(:)=0.0
53
54      ! load high resolution wavenumber data
55      file_id='/stellar_spectra/lam.txt'
56      file_path=TRIM(datadir)//TRIM(file_id)
57
58      open(110,file=file_path,form='formatted',status='old',iostat=ios)
59      if (ios.ne.0) then        ! file not found
60        write(*,*) 'Error from ave_stelspec'
61        write(*,*) 'Data file ',trim(file_id),' not found.'
62        write(*,*)'Check that your path to datagcm:',trim(datadir)
63        write(*,*)' is correct. You can change it in callphys.def with:'
64        write(*,*)' datadir = /absolute/path/to/datagcm'
65        write(*,*)' Also check that there is a ',trim(file_id),' there.'
66        call abort
67      else
68        do ifine=1,Nfine
69          read(110,*) lam(ifine)
70        enddo
71        close(110)
72      endif
73
74      dl=lam(2)-lam(1)
75
76
77      if(stelbbody)then
78         tstellar=stelTbb
79         do ifine=1,Nfine
80            lam_temp=lam(ifine)
81            call blackl(dble(lam_temp*1e-6),dble(tstellar),stel_temp)
82            stel_f(ifine)=stel_temp
83         enddo
84      else
85         ! load high resolution stellar data
86         if(startype.eq.1)then
87            file_id='/stellar_spectra/sol.txt'
88            tstellar=5800.
89         elseif(startype.eq.2)then
90            file_id='/stellar_spectra/gl581.txt'
91            tstellar=3200.
92         elseif(startype.eq.3)then
93            file_id='/stellar_spectra/adleo.txt'
94            tstellar=3200.
95         elseif(startype.eq.4)then
96            file_id='/stellar_spectra/gj644.txt'
97            print*,'Find out tstellar before using this star!'
98            call abort
99         elseif(startype.eq.5)then
100            file_id='/stellar_spectra/hd128167.txt'
101            tstellar=6700. ! Segura et al. (2003)
102         else
103            print*,'Error: unknown star type chosen'
104            call abort
105         endif
106         
107         file_path=TRIM(datadir)//TRIM(file_id)
108         open(111,file=file_path,form='formatted',status='old',iostat=ios)
109         if (ios.ne.0) then        ! file not found
110           write(*,*) 'Error from ave_stelspec'
111           write(*,*) 'Data file ',trim(file_id),' not found.'
112           write(*,*)'Check that your path to datagcm:',trim(datadir)
113           write(*,*)' is correct. You can change it in callphys.def with:'
114           write(*,*)' datadir = /absolute/path/to/datagcm'
115           write(*,*)' Also check that there is a ',trim(file_id),' there.'
116           call abort
117         else
118           do ifine=1,Nfine
119             read(111,*) stel_f(ifine)
120           enddo
121           close(111)
122         endif
123      endif
124
125
126      ! sum data by band
127      band=1
128      do ifine = 1,Nfine
129
130         if(lam(Nfine-ifine+1) .lt. real(10000.0/BWNV(band+1)))then
131            band=band+1
132         endif
133         if(band .gt. L_NSPECTV)then
134            goto 9999 ! ok, ok, I know they're evil
135         endif
136         STELLAR(band)=STELLAR(band)+stel_f(Nfine-ifine+1)*dl
137
138      end do
139
1409999  continue
141      STELLAR=STELLAR/sum(STELLAR)
142
143
144      end subroutine ave_stelspec
Note: See TracBrowser for help on using the repository browser.