source: codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/config/ppsrc/phys/ave_stelspec.f90 @ 224

Last change on this file since 224 was 224, checked in by ymipsl, 10 years ago
File size: 10.1 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!     Generalized to very late spectral types (and Brown dwarfs) Jeremy Leconte (2012)
14!
15!     Called by
16!     ---------
17!     setspv.F
18!     
19!     Calls
20!     -----
21!     none
22!     
23!==================================================================
24
25      use radinc_h, only: L_NSPECTV
26      use radcommon_h, only: BWNV, DWNV, tstellar
27      use datafile_mod, only: datadir
28
29      implicit none
30
31!
32! For Fortran 77/Fortran 90 compliance always use line continuation
33! symbols '&' in columns 73 and 6
34!
35! Group commons according to their type for minimal performance impact
36
37      COMMON/callkeys_l/callrad,corrk,calldifv,UseTurbDiff,calladj      &
38     &   , co2cond,callsoil                                             &
39     &   , season,diurnal,tlocked,rings_shadow,lwrite                   &
40     &   , callstats,calleofdump                                        &
41     &   , enertest                                                     &
42     &   , callgasvis,continuum,H2Ocont_simple,graybody                 &
43     &   , radfixed                                                     &
44     &   , meanOLR, specOLR                                             &
45     &   , kastprof                                                     &
46     &   , nosurf, oblate                                               &     
47     &   , newtonian, testradtimes                                      &
48     &   , check_cpp_match, force_cpp                                   &
49     &   , rayleigh                                                     &
50     &   , stelbbody                                                    &
51     &   , nearco2cond                                                  &
52     &   , tracer, mass_redistrib, varactive, varfixed                  &
53     &   , sedimentation,water,watercond,waterrain                      &
54     &   , aeroco2,aeroh2o,aeroh2so4,aeroback2lay                       &
55     &   , aerofixco2,aerofixh2o                                        &
56     &   , hydrology, sourceevol                                        &
57     &   , CLFvarying                                                   &
58     &   , strictboundcorrk                                             &                                       
59     &   , ok_slab_ocean                                                &
60     &   , ok_slab_sic                                                  &
61     &   , ok_slab_heat_transp                                         
62
63
64      COMMON/callkeys_i/iaervar,iddist,iradia,startype
65     
66      COMMON/callkeys_r/topdustref,Nmix_co2,dusttau,Fat1AU,stelTbb,     &
67     &                  Tstrat,tplanet,obs_tau_col_tropo,               &
68     &                  obs_tau_col_strato,pres_bottom_tropo,           &
69     &                  pres_top_tropo,pres_bottom_strato,              &
70     &                  pres_top_strato,size_tropo,size_strato,satval,  &
71     &                  CLFfixval,n2mixratio,co2supsat,pceil,albedosnow,&
72     &                  maxicethick,Tsaldiff,tau_relax,cloudlvl,        &
73     &                  icetstep,intheat,flatten,Rmean,J2,MassPlanet
74     
75      logical callrad,corrk,calldifv,UseTurbDiff                        &
76     &   , calladj,co2cond,callsoil                                     &
77     &   , season,diurnal,tlocked,rings_shadow,lwrite                   &
78     &   , callstats,calleofdump                                        &
79     &   , callgasvis,continuum,H2Ocont_simple,graybody                 &
80     &   , strictboundcorrk                                             
81
82      logical enertest
83      logical nonideal
84      logical meanOLR
85      logical specOLR
86      logical kastprof
87      logical newtonian
88      logical check_cpp_match
89      logical force_cpp
90      logical testradtimes
91      logical rayleigh
92      logical stelbbody
93      logical ozone
94      logical nearco2cond
95      logical tracer
96      logical mass_redistrib
97      logical varactive
98      logical varfixed
99      logical radfixed
100      logical sedimentation
101      logical water,watercond,waterrain
102      logical aeroco2,aeroh2o,aeroh2so4,aeroback2lay
103      logical aerofixco2,aerofixh2o
104      logical hydrology
105      logical sourceevol
106      logical CLFvarying
107      logical nosurf
108      logical oblate
109      logical ok_slab_ocean
110      logical ok_slab_sic
111      logical ok_slab_heat_transp
112
113      integer iddist
114      integer iaervar
115      integer iradia
116      integer startype
117
118      real topdustref
119      real Nmix_co2
120      real dusttau
121      real Fat1AU
122      real stelTbb
123      real Tstrat
124      real tplanet
125      real obs_tau_col_tropo
126      real obs_tau_col_strato
127      real pres_bottom_tropo
128      real pres_top_tropo
129      real pres_bottom_strato
130      real pres_top_strato
131      real size_tropo
132      real size_strato
133      real satval
134      real CLFfixval
135      real n2mixratio
136      real co2supsat
137      real pceil
138      real albedosnow
139      real maxicethick
140      real Tsaldiff
141      real tau_relax
142      real cloudlvl
143      real icetstep
144      real intheat
145      real flatten
146      real Rmean
147      real J2
148      real MassPlanet
149
150      real*8 STELLAR(L_NSPECTV)
151!      integer startype
152
153      integer Nfine
154      integer,parameter :: Nfineband=200
155      integer ifine
156
157      real,allocatable :: lam(:),stel_f(:)
158      real band,lamm,lamp
159      real dl
160
161      character(len=100)  :: file_id,file_id_lam
162      character(len=200) :: file_path,file_path_lam
163
164      real lam_temp
165      double precision stel_temp
166     
167      integer :: ios ! file opening/reading status
168
169      STELLAR(:)=0.0
170
171      print*,'enter ave_stellspec'
172      if(stelbbody)then
173         tstellar=stelTbb
174         Nfine=L_NSPECTV*Nfineband
175         do band=1,L_NSPECTV
176            lamm=10000.0/BWNV(band+1)
177            lamp=10000.0/BWNV(band)
178            dl=(lamp-lamm)/(Nfineband)
179            do ifine=1,Nfineband
180               lam_temp=lamm+(lamp-lamm)*(ifine-1.)/(Nfineband)
181               call blackl(dble(lam_temp*1e-6),dble(tstellar),stel_temp)
182               STELLAR(band)=STELLAR(band)+stel_temp*dl
183            enddo           
184         end do
185         STELLAR(1:L_NSPECTV)=STELLAR(1:L_NSPECTV)/sum(STELLAR(1:L_NSPECTV))
186      else
187         ! load high resolution stellar data
188         Select Case(startype)
189           Case(1)
190            file_id='/stellar_spectra/sol.txt'
191            tstellar=5800.
192            file_id_lam='/stellar_spectra/lam.txt'
193            Nfine=5000
194           Case(2)
195            file_id='/stellar_spectra/gl581.txt'
196            tstellar=3200.
197            file_id_lam='/stellar_spectra/lam.txt'
198            Nfine=5000
199           Case(3)
200            file_id='/stellar_spectra/adleo.txt'
201            tstellar=3200.
202            file_id_lam='/stellar_spectra/lam.txt'
203            Nfine=5000
204           Case(4)
205            file_id='/stellar_spectra/gj644.txt'
206            print*,'Find out tstellar before using this star!'
207            call abort
208            file_id_lam='/stellar_spectra/lam.txt'
209            Nfine=5000
210           Case(5)
211            file_id='/stellar_spectra/hd128167.txt'
212            tstellar=6700. ! Segura et al. (2003)
213            file_id_lam='/stellar_spectra/lam.txt'
214            Nfine=5000
215           Case(6)
216            file_id='/stellar_spectra/BD_Teff-1600K.txt'
217            tstellar=1600. 
218            file_id_lam='/stellar_spectra/lamBD.txt'
219            Nfine=5000
220           Case(7)
221            file_id='/stellar_spectra/BD_Teff-1000K.txt'
222            tstellar=1000. 
223            file_id_lam='/stellar_spectra/lamBD.txt'
224            Nfine=5000
225           Case(8)
226            file_id='/stellar_spectra/Flux_K5_Teff4700_logg4.5_Met-0.5_BTsettle.dat'
227            tstellar=4700. 
228            file_id_lam='/stellar_spectra/lambda_K5_Teff4700_logg4.5_Met-0.5_BTsettle.dat'
229            Nfine=3986
230           Case Default
231            print*,'Error: unknown star type chosen'
232            call abort
233         End Select
234
235         allocate(lam(Nfine),stel_f(Nfine))
236
237         file_path_lam=TRIM(datadir)//TRIM(file_id_lam)
238         open(110,file=file_path_lam,form='formatted',status='old',iostat=ios)
239         if (ios.ne.0) then        ! file not found
240           write(*,*) 'Error from ave_stelspec'
241           write(*,*) 'Data file ',trim(file_id_lam),' not found.'
242           write(*,*)'Check that your path to datagcm:',trim(datadir)
243           write(*,*)' is correct. You can change it in callphys.def with:'
244           write(*,*)' datadir = /absolute/path/to/datagcm'
245           write(*,*)' Also check that there is a ',trim(file_id_lam),' there.'
246           call abort
247         else
248           do ifine=1,Nfine
249             read(110,*) lam(ifine)
250           enddo
251           close(110)
252         endif
253
254
255         ! load high resolution wavenumber data
256         file_path=TRIM(datadir)//TRIM(file_id)
257         open(111,file=file_path,form='formatted',status='old',iostat=ios)
258         if (ios.ne.0) then        ! file not found
259           write(*,*) 'Error from ave_stelspec'
260           write(*,*) 'Data file ',trim(file_id),' not found.'
261           write(*,*)'Check that your path to datagcm:',trim(datadir)
262           write(*,*)' is correct. You can change it in callphys.def with:'
263           write(*,*)' datadir = /absolute/path/to/datagcm'
264           write(*,*)' Also check that there is a ',trim(file_id),' there.'
265           call abort
266         else
267           do ifine=1,Nfine
268             read(111,*) stel_f(ifine)
269           enddo
270           close(111)
271         endif
272         
273         ! sum data by band
274         band=1
275         Do while(lam(1).lt. real(10000.0/BWNV(band+1)))
276            if (band.gt.L_NSPECTV-1) exit
277            band=band+1
278         enddo
279         dl=lam(2)-lam(1)
280         STELLAR(band)=STELLAR(band)+stel_f(1)*dl
281         do ifine = 2,Nfine
282            if(lam(ifine) .gt. real(10000.0/BWNV(band)))then
283               band=band-1
284            endif
285            if(band .lt. 1) exit
286            dl=lam(ifine)-lam(ifine-1)
287            STELLAR(band)=STELLAR(band)+stel_f(ifine)*dl
288         end do
289               
290         
291         STELLAR(1:L_NSPECTV)=STELLAR(1:L_NSPECTV)/sum(STELLAR(1:L_NSPECTV))
292         deallocate(lam,stel_f)
293         
294      endif
295
296      end subroutine ave_stelspec
Note: See TracBrowser for help on using the repository browser.