libsim  Versione7.2.3
dballe_test.F03
1 #include "config.h"
2 
3 program dballe_test
4 
5 use dballe_class
7 use log4fortran
8 implicit none
9 
10 type(dbasession) :: session, sessionfrom,sessionto
11 type(dbaconnection) :: connection ,connectionfrom
12 integer :: ier, category
13 CHARACTER(len=512) :: a_name
14 
15 !questa chiamata prende dal launcher il nome univoco
16 call l4f_launcher(a_name,a_name_force="dballe_test")
17 !init di log4fortran
18 ier=l4f_init()
19 !imposta a_name
20 category=l4f_category_get(a_name//".main")
21 
22 ! connect to dsn type DBA
23 connection=dbaconnection(dsn="sqlite:dballe_test.sqlite?wipe=true")
24 session=dbasession(connection,wipe=.true.,write=.true.)!,anaflag="write", dataflag="write", attrflag="write")
25 call write1() ! write etherogeneous ensamble of data with attributes and constant data using macro object
26 call write2() ! write an omogeneous vector of data
27 call write3() ! write an etherogeneous ensamble of data
28 call write4() ! write an etherogeneous ensamble of data and attributes
29 call delete1() ! delete one var from the entire DB
30 call delete2() ! delete one var only where are some defined metadati
31 call delete3() ! delete some attributes from one var only where are some defined metadati
32 
33 # ifndef F2003_FULL_FEATURES
34 !close everythings
35 call session%delete()
36 call connection%delete()
37 #endif
38 
39 ! connect to dsn type BUFR file for write
40 session=dbasession(filename="dballe_test.bufr",wipe=.true.,write=.true.,memdb=.false.)
41 call write0() ! write ana on file
42 call write1() ! write etherogeneous ensamble of data with attributes and constant data using macro object
43 call write2() ! write an omogeneous vector of data
44 call write3() ! write an etherogeneous ensamble of data
45 call write4() ! write an etherogeneous ensamble of data and attributes
46 
47 # ifndef F2003_FULL_FEATURES
48 !close everythings
49 call session%delete()
50 #endif
51 
52 ! connect to dsn type BUFR file for read
53 sessionfrom=dbasession(filename="dballe_test.bufr",memdb=.false.)
54 ! connect to dsn type BUFR file for write
55 sessionto=dbasession(filename="dballe_test_copy1f.bufr",wipe=.true.,write=.true.,memdb=.false.)
56 call copy1f() ! copy data and attributes of everythings to an other file
57 
58 # ifndef F2003_FULL_FEATURES
59 !close everythings
60 call sessionto%delete()
61 call sessionfrom%delete()
62 #endif
63 
64 ! use memdb
65 ! connect to dsn type DBA
66 connectionfrom=dbaconnection(dsn="mem:")
67 sessionfrom=dbasession(connectionfrom,wipe=.true.,write=.true.)!,anaflag="write", dataflag="write", attrflag="write")
68 call sessionfrom%messages_open_input(filename="dballe_test.bufr",mode="r",format="BUFR",simplified=.true.)
69 ! connect to dsn type BUFR file for write
70 sessionto=dbasession(filename="dballe_test_copy1fmem.bufr",wipe=.true.,write=.true.,memdb=.false.)
71 call copy1fmem() ! copy data and attributes of everythings to an other file
72 
73 # ifndef F2003_FULL_FEATURES
74 !close everythings
75 call sessionfrom%delete()
76 call connectionfrom%delete()
77 call sessionto%delete()
78 #endif
79 
80 !chiudo il logger
81 CALL l4f_category_delete(category)
82 ier=l4f_fini()
83 
84 contains
85 
86 subroutine write0()
87 
88 type(dbasession) :: sessionana
89 type(dbaanalist) :: anal
90 type(dbaana) :: ana
91 logical :: status
92 
93 ! connect to dsn type BUFR file for write
94 sessionana=dbasession(filename="dballe_test_ana.bufr",wipe=.true.,write=.true.,memdb=.false.)
95 
96 
97 call anal%append(dbaana(lon=11.d0,lat=45.d0))
98 call anal%append(dbaana(lon=12.d0,lat=45.d0))
99 call anal%append(dbaana(lon=13.d0,lat=45.d0))
100 
101 call anal%display()
102 
103 call anal%rewind()
104 !extrude ana
105 do while (anal%element())
106  !call sessionana%extrude(ana=anal%current())
107  ana=anal%current()
108  call ana%extrude(sessionana)
109  call anal%next()
110 end do
111 
112 # ifndef F2003_FULL_FEATURES
113 call sessionana%delete()
114 #endif
115 
116 status=anal%delete()
117 
118 end subroutine write0
119 
120 
121 subroutine write1()
122 
123 type(dbametaanddata),allocatable :: metaanddata(:)
124 type(dbadcv) :: attrv
125 integer :: i
126 
127 print *,"----------------------------------------------"
128 print *,"--------------- write1 ------------------------"
129 
130 allocate(metaanddata(2)) ! one metadata for data and one for constant data
131 
132 metaanddata(1)%metadata=dbametadata( &
133  level=dbalevel(level1=105, l1=2000) &
134  ,timerange=dbatimerange(timerange=4, p1=3600,p2=7200) &
135  ,ana=dbaana(lon=10.d0,lat=45.d0) &
136  ,network=dbanetwork("generic") &
137  ,datetime=dbadatetime(datetime_new(2014,01,06,18,00)))
138 
139 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
140 ! create an etherogeneous ensamble of data
141 allocate (metaanddata(1)%dataattrv%dataattr(2))
142 
143 ! first data
144 allocate (metaanddata(1)%dataattrv%dataattr(1)%dat,source=dbadatai("B13003",85))
145 
146 ! create an etherogeneous ensamble of attr
147 allocate (attrv%dcv(3))
148 allocate (attrv%dcv(1)%dat,source=dbadatar("*B33192",30.))
149 allocate (attrv%dcv(2)%dat,source=dbadatai("*B33193",50))
150 allocate (attrv%dcv(3)%dat,source=dbadatar("*B33194",70.))
151 !assemble data and attribute
152 metaanddata(1)%dataattrv%dataattr(1)%attrv=attrv
153 
154 ! second data
155 allocate (metaanddata(1)%dataattrv%dataattr(2)%dat,source=dbadatai("B12101",27315))
156 ! create an etherogeneous ensamble of attr
157 deallocate(attrv%dcv)
158 allocate (attrv%dcv(2))
159 allocate (attrv%dcv(1)%dat,source=dbadatar("*B33192",30.))
160 allocate (attrv%dcv(2)%dat,source=dbadatai("*B33193",50))
161 !assemble data and attribute
162 metaanddata(1)%dataattrv%dataattr(2)%attrv=attrv
163 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
164 
165 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
166 ! station costant data
167 ! copy the same metadata setting that here we have constant in time data
168 metaanddata(2)%metadata=metaanddata(1)%metadata%dbacontextana()
169 ! create an etherogeneous ensamble of data
170 allocate (metaanddata(2)%dataattrv%dataattr(2))
171 allocate (metaanddata(2)%dataattrv%dataattr(1)%dat,source=dbadatai("B07030",223))
172 allocate (metaanddata(2)%dataattrv%dataattr(1)%attrv%dcv(0)) ! we do not want attributes
173 allocate (metaanddata(2)%dataattrv%dataattr(2)%dat,source=dbadatac("B01019","My beautifull station"))
174 allocate (metaanddata(2)%dataattrv%dataattr(2)%attrv%dcv(0)) ! we do not want attributes
175 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
176 
177 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
178 ! display and save everythings
179 do i=1,size(metaanddata)
180  call metaanddata(i)%display()
181  !call session%extrude(metaanddata=metaanddata(i))
182  call metaanddata(i)%extrude(session)
183 end do
184 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
185 
186 end subroutine write1
187 
188 
189 subroutine write2()
190 
191 type(dbadatar),allocatable :: data(:)
192 
193 type(dbalevel) :: level
194 type(dbatimerange) :: timerange
195 type(dbaana) :: ana
196 type(dbanetwork) :: network
197 type(dbadatetime) :: datetime
198 integer :: i
199 
200 print *,"----------------------------------------------"
201 print *,"--------------- write2 ------------------------"
202 
203 
204 !clear the dballe session
205 call session%unsetall()
206 
207 ! set and display metadata
208 level=dbalevel(level1=105, l1=2000)
209 call level%display()
210 timerange=dbatimerange(timerange=4, p1=3600,p2=7200)
211 call timerange%display()
212 !ana=dbaana(coord=dbacoord(ilon=1000000,ilat=4500000))
213 ana=dbaana(lon=11.d0,lat=45.d0)
214 call ana%display()
215 network=dbanetwork("generic")
216 call network%display()
217 datetime=dbadatetime(datetime_new(2014,01,06,18,00))
218 call datetime%display()
219 
220 ! can set metadata step by step
221 call session%set(level=level)
222 call session%set(timerange=timerange)
223 call session%set(ana=ana)
224 call session%set(network=network)
225 call session%set(datetime=datetime)
226 
227 ! I can use the reverse vision step by step
228 !call level%dbaset(session)
229 !call timerange%dbaset(session)
230 !call ana%dbaset(session)
231 !call network%dbaset(session)
232 !call datetime%dbaset(session)
233 
234 ! create an omogeneous vector of data
235 allocate (data(2),source=[dbadatar("B12102",265.33),dbadatar("B12101",273.15)])
236 
237 !set and display omogeneous data
238 do i =1,size(data)
239  call data(i)%display()
240  call session%set(data=data(i))
241  !call data(i)%dbaset(session)
242 end do
243 
244 !write it in dsn or file
245 call session%prendilo()
246 !!$session%enq(metadata)
247 
248 !close message if I am writing on file
249 call session%close_message()
250 
251 
252 end subroutine write2
253 
254 
255 subroutine write3()
256 
257 type(dbametadata) :: metadata
258 type(dbadcv) :: datav
259 type(dbalevel) :: level
260 type(dbatimerange) :: timerange
261 type(dbaana) :: ana
262 type(dbanetwork) :: network
263 type(dbadatetime) :: datetime
264 
265 print *,"----------------------------------------------"
266 print *,"--------------- write3 ------------------------"
267 
268 
269 !clear the dballe session
270 call session%unsetall()
271 
272 ! set metadata
273 level=dbalevel(level1=105, l1=2000)
274 timerange=dbatimerange(timerange=4, p1=3600,p2=7200)
275 ana=dbaana(lon=12.d0,lat=45.d0)
276 network=dbanetwork("generic")
277 datetime=dbadatetime(datetime_new(2014,01,06,18,00))
278 
279 !assemble metadata
280 metadata=dbametadata(level=level,timerange=timerange,ana=ana,network=network,datetime=datetime)
281 call metadata%display()
282 
283 ! I can set metadata one shot
284 call session%set(metadata)
285 ! or in the reverse vision
286 !call metadata%dbaset(session)
287 
288 call metadata%display()
289 
290 ! create and display an etherogeneous ensamble of data
291 allocate (datav%dcv(2))
292 allocate (datav%dcv(1)%dat,source=dbadatai("B12102",26312))
293 allocate (datav%dcv(2)%dat,source=dbadatar("B12101",273.15))
294 call datav%display()
295 !set data
296 call session%set(datav=datav)
297 
298 ! or in the reverse vision
299 !call datav%dbaset(session)
300 
301 !write it in dsn
302 call session%prendilo()
303 !!$session%enq(metadata)
304 
305 !close message if I am writing on file
306 call session%close_message()
307 
308 end subroutine write3
309 
310 
311 subroutine write4()
312 
313 type(dbametadata) :: metadata
314 type(dbadataattrv) :: dataattrv
315 type(dbadcv) :: attrv
316 type(dbadcv) :: datav
317 
318 print *,"----------------------------------------------"
319 print *,"--------------- write4 ------------------------"
320 
321 ! clear the dballe session
322 call session%unsetall()
323 
324 ! define metadata
325 metadata=dbametadata( &
326  level=dbalevel(level1=105, l1=2000) &
327  ,timerange=dbatimerange(timerange=4, p1=3600,p2=7200) &
328  ,ana=dbaana(lon=13.d0,lat=45.d0) &
329  ,network=dbanetwork("generic") &
330  ,datetime=dbadatetime(datetime_new(2014,01,06,18,00)))
331 
332 call session%set(metadata)
333 
334 ! create and display an etherogeneous ensamble of data
335 allocate (datav%dcv(2))
336 allocate (datav%dcv(1)%dat,source=dbadatai("B12102",26312))
337 allocate (datav%dcv(2)%dat,source=dbadatar("B12101",273.15))
338 
339 ! create and display an etherogeneous ensamble of attr
340 allocate (attrv%dcv(3))
341 allocate (attrv%dcv(1)%dat,source=dbadatai("*B33192",30))
342 allocate (attrv%dcv(2)%dat,source=dbadatac("*B33193","70"))
343 allocate (attrv%dcv(3)%dat,source=dbadatad("*B33194",50.d0))
344 call attrv%display()
345 
346 ! assemble data and attribute
347 allocate (dataattrv%dataattr(2))
348 ! first with attribute
349 allocate (dataattrv%dataattr(1)%dat,source=datav%dcv(1)%dat)
350 dataattrv%dataattr(1)%attrv=attrv
351 ! second without attribute
352 allocate (dataattrv%dataattr(2)%dat,source=datav%dcv(2)%dat)
353 allocate (dataattrv%dataattr(2)%attrv%dcv(0))
354 
355 call dataattrv%display()
356 
357 ! write data and attribute
358 !call session%extrude(dataattrv=dataattrv)
359 call dataattrv%extrude(session)
360 
361 ! work on constant station data
362 call session%set(metadata%dbacontextana())
363 
364 ! write the same data and attribute as constant station data
365 !call session%extrude(dataattrv=dataattrv)
366 call dataattrv%extrude(session)
367 
368 end subroutine write4
369 
370 
371 subroutine delete1()
372 
373 type(dbafilter) :: filter
374 
375 print *,"----------------------------------------------"
376 print *,"--------------- delete1 ----------------------"
377 
378 filter=dbafilter(var="B12101")
379 call session%set(filter=filter)
380 call session%dissolve()
381 
382 end subroutine delete1
383 
384 
385 subroutine delete2()
386 
387 type(dbametadata),allocatable :: metadata(:)
388 type(dbafilter) :: filter
389 
390 print *,"----------------------------------------------"
391 print *,"--------------- delete2 ----------------------"
392 
393 
394 allocate(metadata(1))
395 
396 metadata(1)=dbametadata( &
397  level=dbalevel(level1=105, l1=2000) &
398  ,timerange=dbatimerange(timerange=4, p1=3600,p2=7200) &
399  ,ana=dbaana(lon=11.d0,lat=45.d0) &
400  ,network=dbanetwork("generic") &
401  ,datetime=dbadatetime(datetime_new(2014,01,06,18,00)))
402 
403 
404 filter=dbafilter(var="B12102")
405 call session%set(filter=filter)
406 call session%dissolve(metadata)
407 
408 deallocate(metadata)
409 
410 end subroutine delete2
411 
412 
413 
414 subroutine delete3()
415 
416 type(dbametadata),allocatable :: metadata(:)
417 type(dbafilter) :: filter
418 
419 print *,"----------------------------------------------"
420 print *,"--------------- delete3 ----------------------"
421 
422 
423 allocate(metadata(1))
424 
425 metadata(1)=dbametadata( &
426  level=dbalevel(level1=105, l1=2000) &
427  ,timerange=dbatimerange(timerange=4, p1=3600,p2=7200) &
428  ,ana=dbaana(lon=13.d0,lat=45.d0) &
429  ,network=dbanetwork("generic") &
430  ,datetime=dbadatetime(datetime_new(2014,01,06,18,00)))
431 
432 filter=dbafilter(var="B12102",starvarlist="*B33194,*B33193")
433 call session%set(filter=filter)
434 call session%dissolveattr(metadata)
435 
436 deallocate(metadata)
437 
438 end subroutine delete3
439 
440 
441 subroutine copy1f()
442 type(dbametaanddata), allocatable:: metaanddatav(:)
443 integer :: i
444 
445 print *,"----------------------------------------------"
446 print *,"--------------- copy1f ----------------------"
447 
448 ! data
449 call sessionfrom%filerewind()
450 
451 call sessionfrom%set(filter=dbafilter())
452 call sessionfrom%ingest(metaanddatav)
453 do while (size(metaanddatav) >0)
454  print*,"read/write data; count: ",sessionfrom%count
455  do i =1,size(metaanddatav)
456  print *, "display metaanddatav index: ", i
457 ! call metaanddatav(i)%display()
458  !call sessionto%extrude(metaanddata=metaanddatav(i))
459  call metaanddatav(i)%extrude(sessionto)
460  end do
461 ! call sessionto%close_message()
462  call sessionfrom%ingest(metaanddatav)
463 end do
464 deallocate (metaanddatav)
465 
466 ! constant data
467 call sessionfrom%filerewind()
468 
469 call sessionfrom%set(filter=dbafilter(contextana=.true.))
470 call sessionfrom%ingest(metaanddatav)
471 do while (size(metaanddatav) >0)
472  print*,"read/write data; count: ",sessionfrom%count
473  do i =1,size(metaanddatav)
474  print *, "display metaanddatav index: ", i
475 ! call metaanddatav(i)%display()
476  !call sessionto%extrude(metaanddata=metaanddatav(i))
477  call metaanddatav(i)%extrude(sessionto)
478  end do
479 ! call sessionto%close_message()
480  call sessionfrom%ingest(metaanddatav)
481 end do
482 
483 end subroutine copy1f
484 
485 
486 subroutine copy1fmem()
487 type(dbametaanddata),allocatable :: metaanddatav(:)
488 
489 print *,"----------------------------------------------"
490 print *,"--------------- copy1fmem ----------------------"
491 
492 
493 do while (sessionfrom%messages_read_next())
494 !!$ print*,"read/write message"
495 !!$ call sessionfrom%set(filter=dbafilter())
496 !!$ call sessionfrom%ingest(metaanddatav)
497 !!$ do i =1,size(metaanddatav)
498 !!$ call metaanddatav(i)%display()
499 !!$ call sessionto%extrude(metaanddata=metaanddatav(i))
500 !!$ end do
501 !!$ call sessionto%prendilo()
502 
503  print *,"contextana"
504  call sessionfrom%set(filter=dbafilter(contextana=.true.))
505  call sessionfrom%ingest(metaanddatav)
506 !!$ do i =1,size(metaanddatav)
507 !!$ call metaanddatav(i)%display()
508 !!$ call sessionto%extrude(metaanddata=metaanddatav(i))
509 !!$ end do
510 !!$ call sessionto%close_message()
511  call sessionfrom%remove_all()
512 end do
513 
514 end subroutine copy1fmem
515 
516 end program dballe_test
Class for expressing an absolute time value.
datetime metadata
Classi per la gestione delle coordinate temporali.
manage session handle
character version for dbadata
summ of all metadata pieces
doubleprecision version for dbadata
filter to apply before ingest data
log4fortran destructor
manage connection handle to a DSN
integer version for dbadata
class for import and export data from e to DB-All.e.
classe per la gestione del logging
Global log4fortran constructor.
timerange metadata
real version for dbadata

Generated with Doxygen.