Changeset 3942 for trunk/LMDZ.GENERIC

Timestamp:

Oct 29, 2025, 11:20:13 PM (8 weeks ago)

Author:

mturbet

Message:

remove startype options from ave_stelspec routine

Location:

trunk/LMDZ.GENERIC

Files:

• changelog.txt (modified) (1 diff)
• libf/phystd/ave_stelspec.F90 (modified) (5 diffs)
• libf/phystd/callkeys_mod.F90 (modified) (1 diff)
• libf/phystd/inifis_mod.F90 (modified) (1 diff)

trunk/LMDZ.GENERIC/changelog.txt

@@ -2117 +2117 @@
 Fix some initializations issues in start2archive: call conf_gcm() instead of
 defrun_new() to ensure some defaults (like "timestart") are correctly set.
+
+== 29/10/2025 == MT
+Clean the ave_stelspec routine to remove the startype options
+Now stellar spectra should be specified with stelspec_file and t_stellar options in callphys.def
+Stellar spectra repository now updated in https://web.lmd.jussieu.fr/~lmdz/planets/generic/datagcm/stellar_spectra/

trunk/LMDZ.GENERIC/libf/phystd/ave_stelspec.F90

@@ -12 +12 @@
 !     Robin Wordsworth (2010).
 !     Generalized to very late spectral types (and Brown dwarfs) Jeremy Leconte (2012)
+!     Modified to account for any stellar spectrum file (Lucas Teinturier and Martin Turbet, 2023-2025)
 !
 !     Called by
@@ -26 +27 @@
       use radcommon_h, only: BWNV, DWNV, tstellar
       use datafile_mod, only: datadir
-      use callkeys_mod, only: stelbbody,stelTbb,startype
+      use callkeys_mod, only: stelbbody,stelTbb
       use ioipsl_getin_p_mod, only: getin_p
 
@@ -32 +33 @@
 
       real*8 STELLAR(L_NSPECTV)
-!      integer startype 
 
       integer Nfine
@@ -68 +68 @@
          end do
          STELLAR(1:L_NSPECTV)=STELLAR(1:L_NSPECTV)/sum(STELLAR(1:L_NSPECTV))
-      else
+      else !stelbbody
+         ! look for a " tstellar= ..." option in def files
+         tstellar = -1. ! default
+         call getin_p("tstellar",tstellar) ! default path
+         if (tstellar.eq.-1.) then
+           write(*,*)'Beware that startype is now deprecated, you should use '
+           write(*,*)'stelspec_file and tstellar to define the input stellar spectrum.'
+           write(*,*)'     '
+           write(*,*)'Error: tstellar (effective stellar temperature) needs to be specified'
+           write(*,*)'in callphys.def: tstellar=...'
+           call abort_physic("ave_stelspec", "tstellar needs to be specified",1)
+         end if
+         
+         write(*,*) "Input stellar temperature is:"
+         write(*,*) "tstellar = ",tstellar
+
          ! load high resolution stellar data
-         Select Case(startype)
-           Case(1)
-            file_id='/stellar_spectra/sol.txt'
-            tstellar=5800.
-            file_id_lam='/stellar_spectra/lam.txt'
-            Nfine=5000
-           Case(2)
-            file_id='/stellar_spectra/gl581.txt'
-            tstellar=3200.
-            file_id_lam='/stellar_spectra/lam.txt'
-            Nfine=5000
-           Case(3)
-            file_id='/stellar_spectra/adleo.txt'
-            tstellar=3200.
-            file_id_lam='/stellar_spectra/lam.txt'
-            Nfine=5000
-           Case(4)
-            file_id='/stellar_spectra/gj644.txt'
-            call abort_physic("GJ644", "Tstellar unknown", 1)
-            file_id_lam='/stellar_spectra/lam.txt'
-            Nfine=5000
-           Case(5)
-            file_id='/stellar_spectra/Wasp43_flux.txt'
-            tstellar=4520.
-            file_id_lam='/stellar_spectra/Wasp43_lambda.txt'
-            Nfine=395562
-           Case(6)
-            file_id='/stellar_spectra/BD_Teff-1600K.txt'
-            tstellar=1600. 
-            file_id_lam='/stellar_spectra/lamBD.txt'
-            Nfine=5000
-           Case(7)
-            file_id='/stellar_spectra/BD_Teff-1000K.txt'
-            tstellar=1000. 
-            file_id_lam='/stellar_spectra/lamBD.txt'
-            Nfine=5000
-           Case(8)
-            file_id='/stellar_spectra/Flux_K5_Teff4700_logg4.5_Met-0.5_BTsettle.dat'
-            tstellar=4700. 
-            file_id_lam='/stellar_spectra/lambda_K5_Teff4700_logg4.5_Met-0.5_BTsettle.dat'
-            Nfine=3986
-           Case(9)
-            file_id='/stellar_spectra/Flux_TRAPPIST1.dat'
-            tstellar=2550. 
-            file_id_lam='/stellar_spectra/lambda_TRAPPIST1.dat'
-            Nfine=5000
-           Case(10)
-            file_id='/stellar_spectra/Flux_Proxima.dat'
-            tstellar=3050. 
-            file_id_lam='/stellar_spectra/lambda_Proxima.dat'
-            Nfine=5000
-           Case(11)
-            ! look for a " stelspec_file= ..." option in def files
-            write(*,*) "Input stellar spectra files for radiative transfer is:"
-            stelspec_file = "sun.txt" ! default
-            call getin_p("stelspec_file",stelspec_file) ! default path
-            write(*,*) " stelspec_file = ",trim(stelspec_file)
-            write(*,*) 'Please use ',1,' and only ',1,' header line in ',trim(stelspec_file)
-            ! look for a " tstellar= ..." option in def files
-            write(*,*) "Input stellar temperature is:"
-            tstellar = -1. ! default
-            call getin_p("tstellar",tstellar) ! default path
-            write(*,*) " tstellar = ",tstellar
-            if (tstellar.eq.-1.) then
-              write(*,*)'Error: with startype 11 tstellar need to be specified'
-              write(*,*)'  Specified in callphys.def tstellar=...'
-              call abort_physic("ave_stelspec", "startype 11 tstellar needs to be specified",1)
-            end if
+         ! look for a " stelspec_file= ..." option in def files
+         stelspec_file = "None" ! default
+         call getin_p("stelspec_file",stelspec_file) ! default path
+         
+         write(*,*) "Input stellar spectrum file is:"
+         write(*,*) "stelspec_file = ",trim(stelspec_file)
+         write(*,*) 'Please use ',1,' and only ',1,' header line in ',trim(stelspec_file)
+
+         ! Opening file
+         file_path = trim(datadir)//'/stellar_spectra/'//stelspec_file
+         print*, 'stellar flux : ', file_path
+         OPEN(UNIT=110,FILE=file_path,STATUS='old',iostat=ios)
    
-            ! Opening file
-            file_path = trim(datadir)//'/stellar_spectra/'//stelspec_file
-            print*, 'stellar flux : ', file_path
-            OPEN(UNIT=110,FILE=file_path,STATUS='old',iostat=ios)
-   
-            if (ios /= 0) THEN
-              write(*,*)'Error: cannot open stelspec_file file ', trim(stelspec_file)
-              write(*,*)'It should be in :',trim(datadir),'/stellar_spectra/'
-              write(*,*)'1) You can change the data directory in callphys.def'
-              write(*,*)'   with:'
-              write(*,*)'   datadir=/path/to/the/directory'
-              write(*,*)'2) You can change the input stelspec_file file name in'
-              write(*,*)'   callphys.def with:'
-              write(*,*)'   stelspec_file=filename'
-              call abort_physic("ave_stelspec", "Unable to read stellar flux file", 1)
-            end if
+         if (ios /= 0) THEN
+           write(*,*)'Beware that startype is now deprecated, you should use '
+           write(*,*)'stelspec_file and tstellar to define the input stellar spectrum.'
+           write(*,*)'     '
+           write(*,*)'Error: cannot open stelspec_file file ', trim(stelspec_file)
+           write(*,*)'It should be in :',trim(datadir),'/stellar_spectra/'
+           write(*,*)'1) You can change the data directory in callphys.def'
+           write(*,*)'   with:'
+           write(*,*)'   datadir=/path/to/the/directory'
+           write(*,*)'2) You can change the input stelspec_file file name in'
+           write(*,*)'   callphys.def with:'
+           write(*,*)'   stelspec_file=filename'
+           write(*,*)'You can check the online repository to search for '
+           write(*,*)'available stellar spectra here : '
+           write(*,*)'https://web.lmd.jussieu.fr/~lmdz/planets/generic/datagcm/stellar_spectra/'
+           call abort_physic("ave_stelspec", "Unable to read stellar flux file", 1)
+         end if
 
-            ! Get number of line in the file
-            READ(110,*) ! skip header
-            Nfine = 0
-            do
-              read(110,*,iostat=ios)
-              if (ios<0) exit
-              Nfine = Nfine + 1
-            end do
-            rewind(110) ! Rewind file after counting lines
-            READ(110,*) ! skip header
-           Case Default
-            call abort_physic("ave_stelspec", "Unknown star type chosen", 1)
-         End Select
+         ! Get number of line in the file
+         READ(110,*) ! skip first line header just in case
+         Nfine = 0
+         do
+           read(110,*,iostat=ios)
+           if (ios<0) exit
+           Nfine = Nfine + 1
+         end do
+         rewind(110) ! Rewind file after counting lines
+         READ(110,*) ! skip first line header just in case
 
 !$OMP MASTER
          allocate(lam(Nfine),stel_f(Nfine))
 
-      if (startype.eq.11) then
          do ifine=1,Nfine
            read(110,*) lam(ifine), stel_f(ifine) ! lam [um] stel_f [per unit of wavelength] (integrated and normalized by Fat1AU)
          enddo
-         close(110)
-      else
-         file_path_lam=TRIM(datadir)//TRIM(file_id_lam)
-         open(110,file=file_path_lam,form='formatted',status='old',iostat=ios)
-         if (ios.ne.0) then        ! file not found
-           write(*,*) 'Error from ave_stelspec'
-           write(*,*) 'Data file ',trim(file_id_lam),' not found.'
-           write(*,*)'Check that your path to datagcm:',trim(datadir)
-           write(*,*)' is correct. You can change it in callphys.def with:'
-           write(*,*)' datadir = /absolute/path/to/datagcm'
-           write(*,*)' Also check that there is a ',trim(file_id_lam),' there.'
-           call abort_physic("ave_stelspec", "Unable to find file", 1)
-         else
-           do ifine=1,Nfine
-             read(110,*) lam(ifine)
-           enddo
-           close(110)
-         endif
 
-
-         ! load high resolution wavenumber data
-         file_path=TRIM(datadir)//TRIM(file_id)
-         open(111,file=file_path,form='formatted',status='old',iostat=ios)
-         if (ios.ne.0) then        ! file not found
-           write(*,*) 'Error from ave_stelspec'
-           write(*,*) 'Data file ',trim(file_id),' not found.'
-           write(*,*)'Check that your path to datagcm:',trim(datadir)
-           write(*,*)' is correct. You can change it in callphys.def with:'
-           write(*,*)' datadir = /absolute/path/to/datagcm'
-           write(*,*)' Also check that there is a ',trim(file_id),' there.'
-           call abort_physic("ave_stelspec", "Unable to find file", 1)
-         else
-           do ifine=1,Nfine
-             read(111,*) stel_f(ifine)
-           enddo
-           close(111)
-         endif
-      end if
 !$OMP END MASTER
 !$OMP BARRIER
@@ -243 +162 @@
 !$OMP END MASTER
 !$OMP BARRIER         
-      endif
+      endif !stelbbody
 
       end subroutine ave_stelspec

trunk/LMDZ.GENERIC/libf/phystd/callkeys_mod.F90

@@ -101 +101 @@
       integer,save :: iaervar
       integer,save :: iradia
-      integer,save :: startype
       integer,save :: versH2H2cia
       character(64),save :: H2orthopara_mixture
       integer,save :: nlayaero
-!$OMP THREADPRIVATE(iddist,iaervar,iradia,startype,versH2H2cia,H2orthopara_mixture,nlayaero)
+!$OMP THREADPRIVATE(iddist,iaervar,iradia,versH2H2cia,H2orthopara_mixture,nlayaero)
       integer,dimension(:),allocatable,save :: aeronlay_choice
 !$OMP THREADPRIVATE(aeronlay_choice)

trunk/LMDZ.GENERIC/libf/phystd/inifis_mod.F90

@@ -684 +684 @@
      if (is_master) write(*,*)trim(rname)//": tplanet = ",tplanet
 
-     if (is_master) write(*,*)trim(rname)//": Which star?"
-     startype=1 ! default value = Sol
-     call getin_p("startype",startype)
-     if (is_master) write(*,*)trim(rname)//": startype = ",startype
-
      if (is_master) write(*,*)trim(rname)//": Value of stellar flux at 1 AU?"
      Fat1AU=1356.0 ! default value = Sol today