source: LMDZ5/branches/testing/libf/phylmd/Dust/satellite_out_spla.F90 @ 2791

Last change on this file since 2791 was 2641, checked in by Laurent Fairhead, 8 years ago

Merged trunk changes r2593:2640 into testing branch

File size: 6.1 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
117  INTEGER :: annee, mois, jour
118  REAL :: heure                      !--heure en jour
119  REAL :: pdtphys                    !--pas de temps en seconde
120  REAL :: rlon(1:klon), rlat(1:klon) !--longitude et latitude
121  INTEGER :: masque(1:klon)
122
123  REAL :: J0              !--origine des temps pour les orbites ADEOS
124  PARAMETER (J0=183.91267)
125  REAL :: secinday
126  PARAMETER (secinday=86400.)
127  REAL :: duree_orb       !--Duree d une orbite ADEOS en jour
128  PARAMETER (duree_orb=6055.3715/secinday)
129  REAL :: deltalon        !--Decalage en longitude entre 2 orbites successives
130  PARAMETER (deltalon=41./585.*360.)
131  REAL :: demi_larg_eq    !--demi-largeur d'une orbite a l equateur
132  PARAMETER (demi_larg_eq=11.)
133  REAL :: incli, inclideg
134  PARAMETER (inclideg=98.59)
135  REAL :: RADEG, DTOR, RPI
136  REAL :: demi_periode
137  INTEGER :: jacum(1:12)
138  DATA jacum/0,31,59,90,120,151,181,212,243,273,304,334/
139  INTEGER :: an, orb, i, j
140  REAL :: timepolder, lon0, posnorm, lim_nord, lim_sud
141  REAL :: tempo, lon_cen, demi_larg
142  REAL :: lat_debut, lat_fin, lon_west, lon_east
143  REAL :: zlon            !--rlon mais remis entre 0 et 360
144  REAL :: deltat          !--plage de temps a considerer en jours 
145
146
147  RPI = 4 * atan (1.0)
148  deltat=pdtphys/secinday
149  demi_periode=deltat/2./duree_orb
150  RADEG=180./RPI
151  DTOR=RPI/180.
152  incli=inclideg*DTOR
153
154  an = MOD(annee, 100)
155  timepolder=FLOAT((an-96)*365+jacum(mois)+jour)+heure
156  orb=INT((timepolder-J0)/duree_orb+0.5)
157  lon0=360.-MOD(168.02+FLOAT(orb)*deltalon,360.)
158  posnorm=(timepolder-j0)/duree_orb-FLOAT(orb)
159  j=jacum(mois)+jour
160  lim_nord=60.5 + 25.*(1.-COS((FLOAT(J)+10.)/365.*2.*RPI))
161  lim_sud=-53.0 - 20.*(1.+COS(FLOAT(J)/365.*2*RPI))
162  !--lat de debut
163  lat_debut=MIN(MIN(90.,(-posnorm+demi_periode)*360.),lim_nord)
164  !--lat de fin
165  lat_fin  =MAX(MAX(-90.,(-posnorm-demi_periode)*360.),lim_sud)
166
167  DO i=1, klon
168     masque(i)=0
169     tempo=ASIN( MAX(-1.,MIN(1., -SIN(rlat(i)*DTOR)/SIN(incli))))
170     lon_cen=(ATAN(TAN(tempo)*COS(incli))-duree_orb*tempo)*RADEG
171     demi_larg=demi_larg_eq/COS(rlat(i)*DTOR)
172     IF (ABS(SIN(rlat(i)*DTOR)/SIN(incli)).GE.1.0) demi_larg=200.0
173     IF (rlat(i).GE.lat_fin.AND.rlat(i).LE.lat_debut) THEN
174        IF (demi_larg.GE. 180.) THEN
175           masque(i)=1
176        ELSE 
177           lon_west = MOD(lon0+lon_cen-demi_larg+720., 360.)
178           lon_east = MOD(lon0+lon_cen+demi_larg,      360.)
179           zlon     = MOD(rlon(i)+360.,                360.)
180           IF (lon_west.LE.lon_east) THEN
181              IF (zlon.GE.lon_west.AND.zlon.LE.lon_east) masque(i)=1
182           ELSE               
183              IF (zlon.GE.lon_west.OR.zlon.LE.lon_east) masque(i)=1
184           ENDIF
185        ENDIF
186     ENDIF
187  ENDDO
188
189  RETURN
190END SUBROUTINE swathpolder
191
192
Note: See TracBrowser for help on using the repository browser.