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

Generated with Doxygen.