source: LMDZ5/branches/LMDZ5_SPLA/libf/phylmd/satellite_out_spla.F90

Last change on this file was 2175, checked in by jescribano, 10 years ago

SPLA code included for first time

File size: 6.2 KB
Line 
1
2SUBROUTINE satellite_out_spla(jD_cur,jH_cur,pdtphys,rlat,rlon, &
3                              masque_aqua, masque_terra )
4
5  USE dimphy
6  USE IOIPSL
7  USE write_field_phy
8
9  IMPLICIT NONE
10
11  REAL :: pdtphys, jD_cur,jH_cur, hour
12  INTEGER :: year_cur, mth_cur, day_cur
13  INTEGER :: masque_polder(klon)  ! masque polder
14  INTEGER :: masque_aqua(klon)  ! masque polder
15  INTEGER :: masque_terra(klon)  ! masque polder
16  INTEGER :: i
17  REAL :: overpassaqua, overpassterra
18  REAL,dimension(klon) :: rlat,rlon
19
20
21  masque_polder(:) = 0.
22
23  CALL ju2ymds(jD_cur+jH_cur, year_cur, mth_cur, day_cur, hour)
24!  print*,'JDcur=',jD_cur,'JHcur=',jH_cur,'year_cur' ,year_cur,'mth_cur' ,mth_cur, 'day_cur',day_cur,'hour' ,hour
25
26!  IF ( (year_cur*100.+mth_cur .GE. 199611 ) .AND. (year_cur*100.+mth_cur .LE. 199706)) THEN
27!          CALL swathpolder(year_cur,mth_cur,day_cur,hour/86400., &
28!            pdtphys,rlon,rlat,masque_polder)
29!  ENDIF
30!
31!  DO i=1,klon
32!    IF ( masque_polder(i) .EQ. 1 ) THEN
33!        print *,'polder output point, lon:', rlon(i),', lat: ',rlat(i)
34!    ENDIF
35!  ENDDO
36!  CALL writefield_phy("masque_polder",float(masque_polder),1)
37
38! Aqua
39  masque_aqua(:) = 0.
40  overpassaqua=48600.  ! 13.30 p.m. local time
41  CALL swathpolarsat(year_cur,mth_cur,day_cur,hour, &
42            pdtphys,rlon,rlat,overpassaqua,masque_aqua)
43
44!  DO i=1,klon
45!    IF ( masque_aqua(i) .EQ. 1 ) THEN
46!        print *,'aqua output point, lon:', rlon(i),', lat: ',rlat(i)
47!    ENDIF
48!  ENDDO
49!  CALL writefield_phy("masque_aqua",float(masque_aqua),1)
50
51  masque_terra(:) = 0.
52  overpassterra=37800.  ! 10.30 a.m. local time
53  CALL swathpolarsat(year_cur,mth_cur,day_cur,hour, &
54            pdtphys,rlon,rlat,overpassterra,masque_terra)
55
56!  DO i=1,klon
57!    IF ( masque_terra(i) .EQ. 1 ) THEN
58!        print *,'terra output point, lon:', rlon(i),', lat: ',rlat(i)
59!    ENDIF
60!  ENDDO
61!  CALL writefield_phy("masque_terra",float(masque_terra),1)
62
63  RETURN
64END SUBROUTINE satellite_out_spla
65
66
67
68
69SUBROUTINE swathpolarsat(annee,mois,jour,heure,pdtphys, &
70     rlon,rlat,overpasstime,masque)
71  ! Adaptation from the simple satellite simulator of AeroCom working group Indirect
72  ! forcing ( Johannes Quaas, MPI for Meteorology, Hamburg )
73  ! http://wiki.esipfed.org/index.php/Indirect_forcing
74 
75  USE dimphy
76  IMPLICIT NONE
77
78  INTEGER :: annee, mois, jour, i
79  REAL :: heure                      !--heure en jour
80  REAL :: pdtphys                    !--pas de temps en seconde
81  REAL :: rlon(1:klon), rlat(1:klon) !--longitude et latitude
82  INTEGER :: masque(1:klon)
83  REAL :: localtime(1:klon)
84  REAL :: overpasstime, utctime
85
86  masque(:) = 0
87  utctime=heure
88  DO i=1,klon
89        localtime(i) = utctime + 240. * rlon(i) ! for each degree of longitude east,4 min earlier local time
90!        IF ( localtime(i) > 86400. ) THEN ! this is still the previous day
91!                localtime(i) = localtime(i) - 86400.
92!        ENDIF
93        ! Select 10.30 a.m. ± dt/2
94!       IF ( ABS( localtime(i) - overpasstime ) <= pdtphys/2. ) THEN
95        IF ( ABS(MOD(localtime(i)+86400.*100,86400.) - overpasstime ) <= pdtphys/2. ) THEN
96                masque(i) = 1
97        ENDIF
98   ENDDO
99
100
101END SUBROUTINE swathpolarsat
102
103SUBROUTINE swathpolder(annee,mois,jour,heure,pdtphys, &
104     rlon,rlat,masque)
105!  Adapted from INCA
106
107!  USE inca_dim
108  USE dimphy
109  IMPLICIT NONE
110
111  !--Auteurs : Francois-Marie Breon + Olivier Boucher
112
113  !-- adapted to be used in INCA aerosol module Michael Schulz
114
115  ! not needed?
116  ! #include "dimensions.h"
117  ! #include "../phylmd/dimphy.h"
118!# include "YOMCST_I.h"
119
120  INTEGER :: annee, mois, jour
121  REAL :: heure                      !--heure en jour
122  REAL :: pdtphys                    !--pas de temps en seconde
123  REAL :: rlon(1:klon), rlat(1:klon) !--longitude et latitude
124  INTEGER :: masque(1:klon)
125
126  REAL :: J0              !--origine des temps pour les orbites ADEOS
127  PARAMETER (J0=183.91267)
128  REAL :: secinday
129  PARAMETER (secinday=86400.)
130  REAL :: duree_orb       !--Duree d une orbite ADEOS en jour
131  PARAMETER (duree_orb=6055.3715/secinday)
132  REAL :: deltalon        !--Decalage en longitude entre 2 orbites successives
133  PARAMETER (deltalon=41./585.*360.)
134  REAL :: demi_larg_eq    !--demi-largeur d'une orbite a l equateur
135  PARAMETER (demi_larg_eq=11.)
136  REAL :: incli, inclideg
137  PARAMETER (inclideg=98.59)
138  REAL :: RADEG, DTOR, RPI
139  REAL :: demi_periode
140  INTEGER :: jacum(1:12)
141  DATA jacum/0,31,59,90,120,151,181,212,243,273,304,334/
142  INTEGER :: an, orb, i, j
143  REAL :: timepolder, lon0, posnorm, lim_nord, lim_sud
144  REAL :: tempo, lon_cen, demi_larg
145  REAL :: lat_debut, lat_fin, lon_west, lon_east
146  REAL :: zlon            !--rlon mais remis entre 0 et 360
147  REAL :: deltat          !--plage de temps a considerer en jours 
148
149
150  RPI = 4 * atan (1.0)
151  deltat=pdtphys/secinday
152  demi_periode=deltat/2./duree_orb
153  RADEG=180./RPI
154  DTOR=RPI/180.
155  incli=inclideg*DTOR
156
157  an = MOD(annee, 100)
158  timepolder=FLOAT((an-96)*365+jacum(mois)+jour)+heure
159  orb=INT((timepolder-J0)/duree_orb+0.5)
160  lon0=360.-MOD(168.02+FLOAT(orb)*deltalon,360.)
161  posnorm=(timepolder-j0)/duree_orb-FLOAT(orb)
162  j=jacum(mois)+jour
163  lim_nord=60.5 + 25.*(1.-COS((FLOAT(J)+10.)/365.*2.*RPI))
164  lim_sud=-53.0 - 20.*(1.+COS(FLOAT(J)/365.*2*RPI))
165  !--lat de debut
166  lat_debut=MIN(MIN(90.,(-posnorm+demi_periode)*360.),lim_nord)
167  !--lat de fin
168  lat_fin  =MAX(MAX(-90.,(-posnorm-demi_periode)*360.),lim_sud)
169
170  DO i=1, klon
171     masque(i)=0
172     tempo=ASIN( MAX(-1.,MIN(1., -SIN(rlat(i)*DTOR)/SIN(incli))))
173     lon_cen=(ATAN(TAN(tempo)*COS(incli))-duree_orb*tempo)*RADEG
174     demi_larg=demi_larg_eq/COS(rlat(i)*DTOR)
175     IF (ABS(SIN(rlat(i)*DTOR)/SIN(incli)).GE.1.0) demi_larg=200.0
176     IF (rlat(i).GE.lat_fin.AND.rlat(i).LE.lat_debut) THEN
177        IF (demi_larg.GE. 180.) THEN
178           masque(i)=1
179        ELSE 
180           lon_west = MOD(lon0+lon_cen-demi_larg+720., 360.)
181           lon_east = MOD(lon0+lon_cen+demi_larg,      360.)
182           zlon     = MOD(rlon(i)+360.,                360.)
183           IF (lon_west.LE.lon_east) THEN
184              IF (zlon.GE.lon_west.AND.zlon.LE.lon_east) masque(i)=1
185           ELSE               
186              IF (zlon.GE.lon_west.OR.zlon.LE.lon_east) masque(i)=1
187           ENDIF
188        ENDIF
189     ENDIF
190  ENDDO
191
192  RETURN
193END SUBROUTINE swathpolder
194
195
Note: See TracBrowser for help on using the repository browser.