New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
fbsel.F90 in branches/UKMO/r6232_tracer_advection/NEMOGCM/TOOLS/OBSTOOLS/src – NEMO

source: branches/UKMO/r6232_tracer_advection/NEMOGCM/TOOLS/OBSTOOLS/src/fbsel.F90 @ 9295

Last change on this file since 9295 was 3002, checked in by djlea, 13 years ago

Update documentation for obstools and dataplot. Removal of dataplot code not needed. Addition of headers to some dataplot code. Addition of .exe to command example in obstools.

File size: 8.7 KB
Line 
1PROGRAM fbsel
2   !!---------------------------------------------------------------------
3   !!
4   !!                     ** PROGRAM fbsel **
5   !!
6   !!  ** Purpose : Select or subsample observations
7   !!
8   !!  ** Method  : Use of utilities from obs_fbm.
9   !!
10   !!  ** Action  :
11   !!
12   !!   Usage:
13   !!     fbsel.exe <input filename> <output filename>
14   !!
15   !!   History :
16   !!        ! 2010 (K. Mogensen) Initial version
17   !!----------------------------------------------------------------------
18   USE obs_fbm
19   USE date_utils
20   IMPLICIT NONE
21   TYPE(obfbdata) :: fbdatain,fbdataout
22   CHARACTER(len=256) :: filenamein,filenameout,filenametmp,cnameout
23#ifndef NOIARGCPROTO
24   INTEGER,EXTERNAL :: iargc
25#endif
26   INTEGER,PARAMETER :: maxtyp=1023
27   INTEGER,PARAMETER :: maxdates=20
28   INTEGER :: nqc,ntyp,ndates,ninidate(maxdates),nenddate(maxdates)
29   LOGICAL :: lsplitqc,lsplittyp,lsplitstat
30   INTEGER :: iqc,ityp,idate,istat
31   REAL :: maxlat,minlat,maxlon,minlon
32   CHARACTER(len=ilenwmo) :: cdwmo,cdwmobeg,cdwmoend
33   CHARACTER(len=ilenwmo), DIMENSION(:), POINTER :: clstatids
34   INTEGER :: nstat
35   NAMELIST/namsel/nqc,ntyp,ndates,ninidate,nenddate,lsplitqc,lsplittyp, &
36      &            lsplitstat,maxlat,minlat,maxlon,minlon,cdwmo,&
37      &            cdwmobeg,cdwmoend
38
39   IF (iargc()/=2) THEN
40      WRITE(*,*)'Usage:'
41      WRITE(*,*)'fbsel <input filename> <output filename>'
42      CALL abort
43   ENDIF
44
45   CALL getarg(1,filenamein)
46   CALL getarg(2,filenameout)
47   
48   nqc=-1
49   ntyp=-1
50   ndates=1
51   ninidate=19500101
52   nenddate=21000101
53
54   lsplitqc=.FALSE.
55   lsplittyp=.FALSE.
56   lsplitstat=.FALSE.
57   cdwmo=REPEAT('X',ilenwmo)
58   cdwmobeg=cdwmo
59   cdwmoend=cdwmo
60   maxlat=1e+38
61   minlat=-1e+38
62   maxlon=1e+38
63   minlon=-1e+38
64   OPEN(10,file='namsel.in')
65   READ(10,namsel)
66   CLOSE(10)
67   IF (cdwmobeg==REPEAT('X',ilenwmo)) cdwmobeg=cdwmo
68   IF (cdwmoend==REPEAT('X',ilenwmo)) cdwmoend=cdwmo
69   WRITE(*,namsel)
70
71   CALL init_obfbdata(fbdatain)
72   CALL init_obfbdata(fbdataout)
73
74   WRITE(*,*)'Reading file : ',TRIM(filenamein)
75   CALL read_obfbdata(TRIM(filenamein),fbdatain)
76   WRITE(*,*)'Number of observations before selection = ',fbdatain%nobs
77   DO idate=1,ndates
78      IF (ndates==1) THEN
79         cnameout=filenameout
80      ELSE
81         WRITE(cnameout,'(I2.2,2A)')idate,'_',TRIM(filenameout)
82      ENDIF
83      IF (lsplitqc) THEN
84         DO iqc=1,3
85            CALL fb_sel(fbdatain,fbdataout,iqc,ntyp, &
86               &        ninidate(idate),nenddate(idate), &
87               &        maxlat,minlat,maxlon,minlon,cdwmobeg,cdwmoend)
88            WRITE(filenametmp,'(A,I2.2,A,A)')'qc_',iqc,'_',TRIM(cnameout)
89            IF (fbdataout%nobs>0) THEN
90               WRITE(*,*)'QC selected = ',iqc
91               WRITE(*,*)'Number of observations selected = ',fbdataout%nobs
92               WRITE(*,*)'Writing file : ',TRIM(filenametmp)
93               CALL write_obfbdata(TRIM(filenametmp),fbdataout)
94            ENDIF
95            CALL dealloc_obfbdata(fbdataout)
96         ENDDO
97      ELSEIF (lsplittyp) THEN
98         DO ityp=0,maxtyp
99            CALL fb_sel(fbdatain,fbdataout,nqc,ityp, &
100               &        ninidate(idate),nenddate(idate), &
101               &        maxlat,minlat,maxlon,minlon,cdwmobeg,cdwmoend)
102            WRITE(filenametmp,'(A,I4.4,A,A)')'typ_',ityp,'_',TRIM(cnameout)
103            IF (fbdataout%nobs>0) THEN
104               WRITE(*,*)'Type = ',ityp
105               WRITE(*,*)'Number of observations selected = ',fbdataout%nobs
106               WRITE(*,*)'Writing file : ',TRIM(filenametmp)
107               CALL write_obfbdata(TRIM(filenametmp),fbdataout)
108            ENDIF
109            CALL dealloc_obfbdata(fbdataout)
110         ENDDO
111      ELSEIF (lsplitstat) THEN
112         CALL fb_sel_uniqueids(fbdatain,clstatids,nstat)
113         DO istat=1,nstat
114            CALL fb_sel(fbdatain,fbdataout,nqc,ntyp, &
115               &        ninidate(idate),nenddate(idate), &
116               &        maxlat,minlat,maxlon,minlon,clstatids(istat),clstatids(istat))
117            WRITE(filenametmp,'(4A)')'statid_', &
118               & TRIM(clstatids(istat)),'_',TRIM(cnameout)
119            IF (fbdataout%nobs>0) THEN
120               WRITE(*,*)'Station = ',clstatids(istat)
121               WRITE(*,*)'Number of observations selected = ',fbdataout%nobs
122               WRITE(*,*)'Writing file : ',TRIM(filenametmp)
123               CALL write_obfbdata(TRIM(filenametmp),fbdataout)
124            ENDIF
125            CALL dealloc_obfbdata(fbdataout)
126         ENDDO
127      ELSE
128         CALL fb_sel(fbdatain,fbdataout,nqc,ntyp, & 
129            &        ninidate(idate),nenddate(idate), &
130            &        maxlat,minlat,maxlon,minlon,cdwmobeg,cdwmoend)
131         WRITE(*,*)'Number of observations selected = ',fbdataout%nobs
132         WRITE(*,*)'Writing file : ',TRIM(cnameout)
133         CALL write_obfbdata(TRIM(cnameout),fbdataout)
134         CALL dealloc_obfbdata(fbdataout)
135      ENDIF
136   ENDDO
137
138CONTAINS
139
140   SUBROUTINE fb_sel(fbdatain,fbdataout,nqc,ntyp,ninidate,nenddate,&
141      &              maxlat,minlat,maxlon,minlon,cdwmobeg,cdwmoend)
142      TYPE(obfbdata) :: fbdatain,fbdataout
143      INTEGER :: nqc,ntyp,ninidate,nenddate
144      REAL :: maxlat,minlat,maxlon,minlon
145      CHARACTER(len=ilenwmo) :: cdwmobeg,cdwmoend
146      INTEGER :: jobs
147      INTEGER :: iqc,ityp
148      LOGICAL :: llvalid(fbdatain%nobs)
149      INTEGER :: iyea,imon,iday
150      REAL(KIND=dp) :: zjini,zjend
151      LOGICAL :: lcheckwmo
152
153      lcheckwmo=(cdwmobeg/=REPEAT('X',ilenwmo)).OR.&
154         &      (cdwmoend/=REPEAT('X',ilenwmo))
155      iyea=ninidate/10000
156      imon=ninidate/100-iyea*100
157      iday=ninidate-iyea*10000-imon*100
158      CALL greg2jul(0,0,0,iday,imon,iyea,zjini)
159      iyea=nenddate/10000
160      imon=nenddate/100-iyea*100
161      iday=nenddate-iyea*10000-imon*100
162      CALL greg2jul(0,0,0,iday,imon,iyea,zjend)
163      DO jobs = 1, fbdatain%nobs
164         llvalid(jobs)=.TRUE.
165         IF (nqc/=-1) THEN
166            CALL check_prof(fbdatain,jobs,iqc)
167            llvalid(jobs)=(iqc==nqc).AND.llvalid(jobs)
168         ENDIF
169         IF (ntyp/=-1) THEN
170            READ(fbdatain%cdtyp(jobs),'(I4)')ityp
171            llvalid(jobs)=(ityp==ntyp).AND.llvalid(jobs)
172         ENDIF
173         IF (ninidate/=-1) THEN
174            llvalid(jobs)=(fbdatain%ptim(jobs)>zjini).AND.llvalid(jobs)
175         ENDIF
176         IF (nenddate/=-1) THEN
177            llvalid(jobs)=(fbdatain%ptim(jobs)<=zjend).AND.llvalid(jobs)
178         ENDIF
179         llvalid(jobs)=(fbdatain%pphi(jobs)<=maxlat).AND.       &
180            &          (fbdatain%pphi(jobs)>=minlat).AND.       &
181            &          (((fbdatain%plam(jobs)<=maxlon).AND.     &
182            &            (fbdatain%plam(jobs)>=minlon)).OR.     &
183            &           ((fbdatain%plam(jobs)+360<=maxlon).AND. &
184            &            (fbdatain%plam(jobs)+360>=minlon)).OR. &
185            &           ((fbdatain%plam(jobs)-360<=maxlon).AND. &
186            &            (fbdatain%plam(jobs)-360>=minlon))).AND.llvalid(jobs) 
187         IF (lcheckwmo) THEN
188            llvalid(jobs)=LGE(TRIM(fbdatain%cdwmo(jobs)),TRIM(cdwmobeg)) &
189               &    .AND. LLE(TRIM(fbdatain%cdwmo(jobs)),TRIM(cdwmoend)) &
190               &    .AND. llvalid(jobs)
191         ENDIF
192         ! Add more checks here...
193      ENDDO
194
195      CALL subsamp_obfbdata(fbdatain,fbdataout,llvalid)
196
197   END SUBROUTINE fb_sel
198
199   SUBROUTINE fb_sel_uniqueids(fbdatain,clstatids,nstat)
200      TYPE(obfbdata) :: fbdatain
201      CHARACTER(len=ilenwmo), DIMENSION(:), POINTER :: clstatids
202      INTEGER :: nstat
203      INTEGER :: jobs,kobs
204      LOGICAL, DIMENSION(fbdatain%nobs) :: lunique
205
206      lunique(:)=.TRUE.
207      DO jobs=1,fbdatain%nobs
208         IF (.NOT.lunique(jobs)) CYCLE
209         DO kobs=jobs+1,fbdatain%nobs
210            IF (.NOT.lunique(kobs)) CYCLE
211            IF (fbdatain%cdwmo(jobs)==fbdatain%cdwmo(kobs)) THEN
212               lunique(kobs)=.FALSE.
213            ENDIF
214         ENDDO
215      ENDDO
216      nstat=COUNT(lunique)
217      ALLOCATE(clstatids(nstat))
218      kobs=0
219      DO jobs=1,fbdatain%nobs
220         IF (lunique(jobs)) THEN
221            kobs=kobs+1
222            clstatids(kobs)=fbdatain%cdwmo(jobs)
223         ENDIF
224      ENDDO
225      WRITE(*,*)'Unique station ids'
226      DO jobs=1,nstat
227         WRITE(*,'(I5,1X,A)')jobs,clstatids(jobs)
228      ENDDO
229
230   END SUBROUTINE fb_sel_uniqueids
231   
232   SUBROUTINE check_prof(fbdata,iobs,iqc)
233     
234      TYPE(obfbdata) :: fbdata
235      INTEGER :: iobs,iqc
236      INTEGER :: i,ivar
237     
238      LOGICAL :: lpart,lfull
239
240      lpart=.false.
241      lfull=.true.
242      DO ivar=1,fbdata%nvar
243         DO i=1,fbdata%nlev
244            IF ((fbdata%ivlqc(i,iobs,ivar)>2).AND.&
245               &(fbdata%ivlqc(i,iobs,ivar)<9)) lpart = .TRUE.
246            IF (fbdata%ivlqc(i,iobs,ivar)<=2) lfull = .FALSE.
247         ENDDO
248      ENDDO
249
250      IF(lfull) THEN
251         iqc=3
252      ELSEIF (lpart) then
253         iqc=2
254      ELSE
255         iqc=1
256      ENDIF
257
258   END SUBROUTINE check_prof
259
260END PROGRAM fbsel
Note: See TracBrowser for help on using the repository browser.