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

Last change on this file since 224 was 224, checked in by ymipsl, 10 years ago
File size: 9.2 KB
Line 
1      subroutine setspv
2
3!==================================================================
4!     
5!     Purpose
6!     -------
7!     Set up spectral intervals, stellar spectrum and Rayleigh
8!     opacity in the shortwave.
9!     
10!     Authors
11!     -------
12!     Adapted from setspv in the NASA Ames radiative code by
13!     Robin Wordsworth (2009).
14!
15!     Called by
16!     ---------
17!     callcorrk.F
18!     
19!     Calls
20!     -----
21!     ave_stelspec.F
22!     
23!==================================================================
24
25      use radinc_h,    only: L_NSPECTV, corrkdir, banddir
26      use radcommon_h, only: BWNV,BLAMV,WNOV,DWNV,WAVEV, &
27                             STELLARF,TAURAY
28      use datafile_mod, only: datadir
29
30      implicit none
31
32!-----------------------------------------------------------------------
33! INCLUDE "comcstfi.h"
34
35      common/comcstfi/pi,rad,g,r,cpp,rcp,dtphys,daysec,mugaz,omeg
36      common/comcstfi/avocado!,molrad,visc
37     
38      real pi,rad,g,r,cpp,rcp,dtphys,daysec,mugaz,omeg
39      real avocado!,molrad,visc
40
41!
42! For Fortran 77/Fortran 90 compliance always use line continuation
43! symbols '&' in columns 73 and 6
44!
45! Group commons according to their type for minimal performance impact
46
47      COMMON/callkeys_l/callrad,corrk,calldifv,UseTurbDiff,calladj      &
48     &   , co2cond,callsoil                                             &
49     &   , season,diurnal,tlocked,rings_shadow,lwrite                   &
50     &   , callstats,calleofdump                                        &
51     &   , enertest                                                     &
52     &   , callgasvis,continuum,H2Ocont_simple,graybody                 &
53     &   , radfixed                                                     &
54     &   , meanOLR, specOLR                                             &
55     &   , kastprof                                                     &
56     &   , nosurf, oblate                                               &     
57     &   , newtonian, testradtimes                                      &
58     &   , check_cpp_match, force_cpp                                   &
59     &   , rayleigh                                                     &
60     &   , stelbbody                                                    &
61     &   , nearco2cond                                                  &
62     &   , tracer, mass_redistrib, varactive, varfixed                  &
63     &   , sedimentation,water,watercond,waterrain                      &
64     &   , aeroco2,aeroh2o,aeroh2so4,aeroback2lay                       &
65     &   , aerofixco2,aerofixh2o                                        &
66     &   , hydrology, sourceevol                                        &
67     &   , CLFvarying                                                   &
68     &   , strictboundcorrk                                             &                                       
69     &   , ok_slab_ocean                                                &
70     &   , ok_slab_sic                                                  &
71     &   , ok_slab_heat_transp                                         
72
73
74      COMMON/callkeys_i/iaervar,iddist,iradia,startype
75     
76      COMMON/callkeys_r/topdustref,Nmix_co2,dusttau,Fat1AU,stelTbb,     &
77     &                  Tstrat,tplanet,obs_tau_col_tropo,               &
78     &                  obs_tau_col_strato,pres_bottom_tropo,           &
79     &                  pres_top_tropo,pres_bottom_strato,              &
80     &                  pres_top_strato,size_tropo,size_strato,satval,  &
81     &                  CLFfixval,n2mixratio,co2supsat,pceil,albedosnow,&
82     &                  maxicethick,Tsaldiff,tau_relax,cloudlvl,        &
83     &                  icetstep,intheat,flatten,Rmean,J2,MassPlanet
84     
85      logical callrad,corrk,calldifv,UseTurbDiff                        &
86     &   , calladj,co2cond,callsoil                                     &
87     &   , season,diurnal,tlocked,rings_shadow,lwrite                   &
88     &   , callstats,calleofdump                                        &
89     &   , callgasvis,continuum,H2Ocont_simple,graybody                 &
90     &   , strictboundcorrk                                             
91
92      logical enertest
93      logical nonideal
94      logical meanOLR
95      logical specOLR
96      logical kastprof
97      logical newtonian
98      logical check_cpp_match
99      logical force_cpp
100      logical testradtimes
101      logical rayleigh
102      logical stelbbody
103      logical ozone
104      logical nearco2cond
105      logical tracer
106      logical mass_redistrib
107      logical varactive
108      logical varfixed
109      logical radfixed
110      logical sedimentation
111      logical water,watercond,waterrain
112      logical aeroco2,aeroh2o,aeroh2so4,aeroback2lay
113      logical aerofixco2,aerofixh2o
114      logical hydrology
115      logical sourceevol
116      logical CLFvarying
117      logical nosurf
118      logical oblate
119      logical ok_slab_ocean
120      logical ok_slab_sic
121      logical ok_slab_heat_transp
122
123      integer iddist
124      integer iaervar
125      integer iradia
126      integer startype
127
128      real topdustref
129      real Nmix_co2
130      real dusttau
131      real Fat1AU
132      real stelTbb
133      real Tstrat
134      real tplanet
135      real obs_tau_col_tropo
136      real obs_tau_col_strato
137      real pres_bottom_tropo
138      real pres_top_tropo
139      real pres_bottom_strato
140      real pres_top_strato
141      real size_tropo
142      real size_strato
143      real satval
144      real CLFfixval
145      real n2mixratio
146      real co2supsat
147      real pceil
148      real albedosnow
149      real maxicethick
150      real Tsaldiff
151      real tau_relax
152      real cloudlvl
153      real icetstep
154      real intheat
155      real flatten
156      real Rmean
157      real J2
158      real MassPlanet
159
160      logical file_ok
161
162      integer N, M, file_entries
163
164      character(len=30)  :: temp1
165      character(len=200) :: file_id
166      character(len=200) :: file_path
167
168      real*8 :: lastband(2)
169
170      real*8 STELLAR(L_NSPECTV)
171      real*8 sum, dummy
172
173      !! used to count lines
174      integer :: nb
175      integer :: ierr
176
177!=======================================================================
178!     Set up spectral bands - wavenumber [cm^(-1)]. Go from smaller to
179!     larger wavenumbers, the same as in the IR.
180
181      write(temp1,'(i2.2)') L_NSPECTV
182      file_id='/corrk_data/'//trim(adjustl(banddir))//'/narrowbands_VI.in' 
183      file_path=TRIM(datadir)//TRIM(file_id)
184
185      ! check that the file exists
186      inquire(FILE=file_path,EXIST=file_ok)
187      if(.not.file_ok) then
188         write(*,*)'The file ',TRIM(file_path)
189         write(*,*)'was not found by setspv.F90, exiting.'
190         write(*,*)'Check that your path to datagcm:',trim(datadir)
191         write(*,*)' is correct. You can change it in callphys.def with:'
192         write(*,*)' datadir = /absolute/path/to/datagcm'
193         write(*,*)'Also check that the corrkdir you chose in callphys.def exists.'
194         call abort
195      endif
196   
197      nb=0
198      ierr=0
199      ! check that the file contains the right number of bands
200      open(131,file=file_path,form='formatted')
201      read(131,*,iostat=ierr) file_entries
202      do while (ierr==0)
203        read(131,*,iostat=ierr) dummy
204        if (ierr==0) nb=nb+1
205      enddo
206      close(131)
207
208      write(*,*) 'setspv: L_NSPECTV = ',L_NSPECTV, 'in the model '
209      write(*,*) '        there are   ',nb, 'entries in ',TRIM(file_path)
210      if(nb.ne.L_NSPECTV) then
211         write(*,*) 'MISMATCH !! I stop here'
212         call abort
213      endif
214
215      ! load and display the data
216      open(111,file=file_path,form='formatted')
217      read(111,*) 
218       do M=1,L_NSPECTV-1
219         read(111,*) BWNV(M)
220      end do
221      read(111,*) lastband
222      close(111)
223      BWNV(L_NSPECTV)  =lastband(1)
224      BWNV(L_NSPECTV+1)=lastband(2)
225
226
227      print*,'setspv: VI band limits:'
228      do M=1,L_NSPECTV+1
229         print*,m,'-->',BWNV(M),' cm^-1'
230      end do
231      print*,' '
232
233!     Set up mean wavenumbers and wavenumber deltas.  Units of
234!     wavenumbers is cm^(-1); units of wavelengths is microns.
235
236      do M=1,L_NSPECTV
237         WNOV(M)  = 0.5*(BWNV(M+1)+BWNV(M))
238         DWNV(M)  = BWNV(M+1)-BWNV(M)
239         WAVEV(M) = 1.0E+4/WNOV(M)
240         BLAMV(M) = 0.01/BWNV(M)
241      end do
242      BLAMV(M) = 0.01/BWNV(M) ! wavelength in METERS for aerosol stuff
243!     note M=L_NSPECTV+1 after loop due to Fortran bizarreness
244
245!=======================================================================
246!     Set up stellar spectrum
247
248      write(*,*)'setspv: Interpolating stellar spectrum from the hires data...'
249      call ave_stelspec(STELLAR)
250
251!     Sum the stellar flux, and write out the result. 
252      sum = 0.0 
253      do N=1,L_NSPECTV
254         STELLARF(N) = STELLAR(N) * Fat1AU
255         sum         = sum+STELLARF(N)
256      end do
257      write(6,'("setspv: Stellar flux at 1 AU = ",f7.2," W m-2")') sum
258      print*,' '
259
260
261!=======================================================================
262!     Set up the wavelength independent part of the Rayleigh scattering.
263!     The pressure dependent part will be computed elsewhere (OPTCV).
264!     WAVEV is in microns.  There is no Rayleigh scattering in the IR.
265
266      if(rayleigh) then
267         call calc_rayleigh
268      else
269         print*,'setspv: No Rayleigh scattering, check for NaN in output!'
270         do N=1,L_NSPECTV
271            TAURAY(N) = 1E-16
272         end do
273      endif
274
275      RETURN
276    END subroutine setspv
Note: See TracBrowser for help on using the repository browser.