libsim Versione 7.2.4
example_dballe.F03

Sample program to demostrate the dballe_class module.

Sample program to demostrate the dballe_class module. This module have examples to read/write/manipulate data from/to DB or BUFR.

1program example_dballe
2
6implicit none
7
8integer :: i
9type(dbasession) :: session,sessionfrom,sessionto
10type(dbaconnection) :: connection,connectionfrom
11
12
13type (dbacoord) :: coord
14type (dbaana) :: ana
15type (dbadatetime) :: mydatetime
16type (dbalevel) :: level
17type (dbanetwork) :: network
18type (dbatimerange) :: timerange
19type (dbametadata) :: metadata
20type (dbadatar) :: datar
21type (dbametaanddatar) :: metaanddatar
22type (dbadc) :: dc
23type (dbadcv) :: dcv
24type (dbametaanddatav) :: metaanddatav
25type (dbadataattrv) :: dataattrv
26type (dbametaanddata) :: metaandata
27
28
29integer :: category,ier
30character(len=512):: a_name
31
32!questa chiamata prende dal launcher il nome univoco
33call l4f_launcher(a_name,a_name_force="example",a_name_append="main")
34!init di log4fortran
35ier=l4f_init()
36!imposta a_name
37category=l4f_category_get(a_name//".main")
38
39call l4f_category_log(category,l4f_info,"example_dballe")
40
41!Oggetti:
42
43!dbacoord : fortran 2003 interface to geo_coord
44coord=dbacoord(lon=10.d0,lat=45.d0)
45!dbaana : ana metadata
46ana=dbaana(coord,ident="mianave")
47
48!dbadatetime : datetime metadata
49mydatetime=dbadatetime(datetime_new(2014,01,06,18,00))
50
51!dbalevel : level metadata
52level=dbalevel(level1=105, l1=2000)
53
54!dbanetwork : network metadata
55network=dbanetwork("generic")
56
57!dbatimerange : timerange metadata
58timerange=dbatimerange(timerange=4, p1=3600,p2=7200)
59
60!dbametadata : summ of all metadata pieces
61metadata=dbametadata(level=level,timerange=timerange,ana=ana,network=network,datetime=mydatetime)
62
63!dbametaanddatar : metadata and real data
64metaanddatar%metadata=metadata
65metaanddatar%dbadatar=datar
66!dbametaanddatarlist : metadata and real data double linked list
67
68!dbametaanddatai : metadata and integer data
69!dbametaanddatailist : metadata and integer data double linked list
70
71!dbametaanddatab : metadata and byte data
72!dbametaanddatablist : metadata and byte data double linked list
73
74!dbametaanddatac : metadata and character data
75!dbametaanddataclist : metadata and character data double linked list
76
77!dbametaanddatad : metadata and doubleprecision data
78!dbametaanddatadlist : metadata and diubleprecision data double linked list
79
80!dbadc : container for dbadata (used for promiscuous vector of data)
81allocate(dc%dat,source=dbadatai("B12102",26312))
82
83!dbadcv : vector of container of dbadata
84allocate (dcv%dcv(2))
85allocate (dcv%dcv(1)%dat,source=dbadatai("B12102",26312))
86allocate (dcv%dcv(2)%dat,source=dbadatar("B12101",273.15))
87
88!dbametaanddatav one metadata plus vector of container of dbadata
89metaanddatav%metadata=metadata
90metaanddatav%datav=dcv
91
92!dbadataattrv : vector of dbadataattr (more data plus attributes)
93allocate(dataattrv%dataattr(2))
94dataattrv%dataattr(1)%dbadc=dc
95dataattrv%dataattr(1)%attrv=dcv
96
97!dbametaanddata : one metadata plus vector of container of dbadata plus attributes
98metaandata%metadata=metadata
99metaandata%dataattrv=dataattrv
100
101!dbametaanddatalist double linked list of dbametaanddata
102
103
104call l4f_category_log(category,l4f_info,"connect to dsn type DBA")
105
106connection=dbaconnection(dsn="sqlite:/tmp/dballe.sqlite")
107session=dbasession(connection,wipe=.true.,anaflag="write", dataflag="write", attrflag="write")
108
109call l4f_category_log(category,l4f_info,"write1")
110call write1() ! write etherogeneous ensamble of data with attributes and constant data using macro object
111call l4f_category_log(category,l4f_info,"write2")
112call write2() ! write an omogeneous vector of data
113call l4f_category_log(category,l4f_info,"write3")
114call write3() ! write an etherogeneous ensamble of data
115call l4f_category_log(category,l4f_info,"write4")
116call write4() ! write an etherogeneous ensamble of data and attributes
117call l4f_category_log(category,l4f_info,"write5")
118call write5() ! write a content of DSN to a file with filter
119call l4f_category_log(category,l4f_info,"read2")
120call read2() ! read data and attributes for data filtered and ordered by btable with type predefined
121call l4f_category_log(category,l4f_info,"delete1")
122call delete1() ! delete one var from the entire DB
123call l4f_category_log(category,l4f_info,"delete2")
124call delete2() ! delete one var only where are some defined metadata
125call l4f_category_log(category,l4f_info,"delete3")
126call delete3() ! delete some attributes from one var only where are some defined metadata
127call l4f_category_log(category,l4f_info,"read1")
128call read1() ! read ana, data and attributes for constant data, data and attributes for data
129call l4f_category_log(category,l4f_info,"read3")
130call read3() ! read an omogeneous vector of data
131call l4f_category_log(category,l4f_info,"read4")
132call read4l() ! read data and attributes for data in a double linked list
133call l4f_category_log(category,l4f_info,"copy1")
134call copy1() ! copy data and attributes of everythings to an other dsn
135call l4f_category_log(category,l4f_info,"copy2")
136call copy2() ! copy data and constant data to an other dsn
137
138!close everythings
139call session%delete()
140call connection%delete()
141
142
143call l4f_category_log(category,l4f_info,"connect to dsn type BUFR file for write")
144session=dbasession(filename="example_dballe.bufr",wipe=.true.,write=.true.,memdb=.false.)
145
146call l4f_category_log(category,l4f_info,"write0")
147call write0() ! write ana on file
148call l4f_category_log(category,l4f_info,"write1")
149call write1() ! write etherogeneous ensamble of data with attributes and constant data using macro object
150call l4f_category_log(category,l4f_info,"write2")
151call write2() ! write an omogeneous vector of data
152call l4f_category_log(category,l4f_info,"write3")
153call write3() ! write an etherogeneous ensamble of data
154call l4f_category_log(category,l4f_info,"write4")
155call write4() ! write an etherogeneous ensamble of data and attributes
156
157!close everythings
158call session%delete()
159
160
161call l4f_category_log(category,l4f_info,"connect to dsn type BUFR file for read")
162session=dbasession(filename="example_dballe.bufr",memdb=.false.)
163
164call l4f_category_log(category,l4f_info,"read1lf")
165call read1lf() ! read ana, data and attributes for constant data, data and attributes for data in list
166call l4f_category_log(category,l4f_info,"read2lf")
167call read2lf() ! read data and attributes for data filtered and ordered by btable with type predefined in list
168
169!! note: I cannot read from file without filter (filter do not work on file) and put everythings in real matrix
170!! we get error putting somethings that do not fit in real (like station name)
171!! call read3lf() ! read an omogeneous vector of data
172
173! use memdb
174call l4f_category_log(category,l4f_info,"read3lfmem")
175call read3lfmem() ! read an omogeneous vector of data from file using list and mem (and filters)
176call l4f_category_log(category,l4f_info,"readmem")
177call readmem() ! read an list of data from file using mem (and filters)
178
179call l4f_category_log(category,l4f_info,"connect to dsn type BUFR file for read")
180sessionfrom=dbasession(filename="example_dballe.bufr",memdb=.false.)
181call l4f_category_log(category,l4f_info,"connect to dsn type BUFR file for write")
182sessionto=dbasession(filename="example_dballe_copy1f.bufr",wipe=.true.,write=.true.,memdb=.false.)
183call l4f_category_log(category,l4f_info,"copy1f")
184call copy1f() ! copy data and attributes of everythings to an other file
185! todo
186!call l4f_category_log(category,L4F_INFO,"copy2f")
187!call copy2f() ! copy data and constant data to an other file
188!close everythings
189call sessionto%delete()
190call sessionfrom%delete()
191
192
193call l4f_category_log(category,l4f_info,"use memdb: connect to dsn type DBA")
194connectionfrom=dbaconnection(dsn="mem:")
195sessionfrom=dbasession(connectionfrom,wipe=.true.,anaflag="write", dataflag="write", attrflag="write")
196call sessionfrom%messages_open_input(filename="example_dballe.bufr",mode="r",format="BUFR",simplified=.true.)
197
198call l4f_category_log(category,l4f_info,"connect to dsn type BUFR file for write")
199sessionto=dbasession(filename="example_dballe_copy1fmem.bufr",wipe=.true.,write=.true.,memdb=.false.)
200call l4f_category_log(category,l4f_info,"copy1fmem")
201call copy1fmem() ! copy data and attributes of everythings to an other file using mem for input
202call sessionto%delete()
203
204call sessionfrom%messages_open_output(filename="example_dballe_copyf2mem.bufr")
205call l4f_category_log(category,l4f_info,"copy1f2mem")
206call copy1f2mem() ! copy data and attributes of everythings to an other file using mem for input and output
207
208!close everythings
209call sessionfrom%delete()
210call connectionfrom%delete()
211
212contains
213
214subroutine write0()
215
216type(dbasession) :: sessionana
217type(dbaanalist) :: anal
218type(dbaana) :: ana
219logical :: status
220
221! connect to dsn type BUFR file for write
222sessionana=dbasession(filename="example_dballe_ana.bufr",wipe=.true.,write=.true.,memdb=.false.)
223
224
225call anal%append(dbaana(lon=11.d0,lat=45.d0))
226call anal%append(dbaana(lon=12.d0,lat=45.d0))
227call anal%append(dbaana(lon=13.d0,lat=45.d0))
228
229call anal%display()
230
231call anal%rewind()
232!extrude ana
233do while (anal%element())
234 ana=anal%current()
235 call ana%extrude(sessionana)
236 call anal%next()
237end do
238
239call sessionana%delete()
240status=anal%delete()
241
242end subroutine write0
243
244
245subroutine write1()
246
247type(dbametaanddata),allocatable :: metaanddata(:)
248type(dbadcv) :: attrv
249
250print *,"----------------------------------------------"
251print *,"--------------- write1 ------------------------"
252
253allocate(metaanddata(2)) ! one metadata for data and one for constant data
254
255metaanddata(1)%metadata=dbametadata( &
256 level=dbalevel(level1=105, l1=2000) &
257 ,timerange=dbatimerange(timerange=4, p1=3600,p2=7200) &
258 ,ana=dbaana(lon=10.d0,lat=45.d0) &
259 ,network=dbanetwork("generic") &
260 ,datetime=dbadatetime(datetime_new(2014,01,06,18,00)))
261
262!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
263! create an etherogeneous ensamble of data
264allocate (metaanddata(1)%dataattrv%dataattr(2))
265
266! first data
267allocate (metaanddata(1)%dataattrv%dataattr(1)%dat,source=dbadatai("B13003",85))
268
269! create an etherogeneous ensamble of attr
270allocate (attrv%dcv(3))
271allocate (attrv%dcv(1)%dat,source=dbadatar("*B33192",30.))
272allocate (attrv%dcv(2)%dat,source=dbadatai("*B33193",50))
273allocate (attrv%dcv(3)%dat,source=dbadatar("*B33194",70.))
274!assemble data and attribute
275metaanddata(1)%dataattrv%dataattr(1)%attrv=attrv
276
277! second data
278allocate (metaanddata(1)%dataattrv%dataattr(2)%dat,source=dbadatai("B12101",27315))
279! create an etherogeneous ensamble of attr
280deallocate(attrv%dcv)
281allocate (attrv%dcv(2))
282allocate (attrv%dcv(1)%dat,source=dbadatar("*B33192",30.))
283allocate (attrv%dcv(2)%dat,source=dbadatai("*B33193",50))
284!assemble data and attribute
285metaanddata(1)%dataattrv%dataattr(2)%attrv=attrv
286!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
287
288!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
289! station costant data
290! copy the same metadata setting that here we have constant in time data
291metaanddata(2)%metadata=metaanddata(1)%metadata%dbacontextana()
292! create an etherogeneous ensamble of data
293allocate (metaanddata(2)%dataattrv%dataattr(2))
294allocate (metaanddata(2)%dataattrv%dataattr(1)%dat,source=dbadatai("B07030",223))
295allocate (metaanddata(2)%dataattrv%dataattr(1)%attrv%dcv(0)) ! we do not want attributes
296allocate (metaanddata(2)%dataattrv%dataattr(2)%dat,source=dbadatac("B01019","My beautifull station"))
297allocate (metaanddata(2)%dataattrv%dataattr(2)%attrv%dcv(0)) ! we do not want attributes
298!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
299
300!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
301! display and save everythings
302do i=1,size(metaanddata)
303 call metaanddata(i)%display()
304 !call session%extrude(metaanddata(i))
305 call metaanddata(i)%extrude(session)
306end do
307!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
308
309end subroutine write1
310
311
312subroutine write2()
313
314type(dbadatar),allocatable :: data(:)
315
316type(dbalevel) :: level
317type(dbatimerange) :: timerange
318type(dbaana) :: ana
319type(dbanetwork) :: network
320type(dbadatetime) :: datetime
321
322print *,"----------------------------------------------"
323print *,"--------------- write2 ------------------------"
324
325
326!clear the dballe session
327call session%unsetall()
328
329! set and display metadata
330level=dbalevel(level1=105, l1=2000)
331call level%display()
332timerange=dbatimerange(timerange=4, p1=3600,p2=7200)
333call timerange%display()
334!ana=dbaana(coord=dbacoord(ilon=1000000,ilat=4500000))
335ana=dbaana(lon=11.d0,lat=45.d0)
336call ana%display()
337network=dbanetwork("generic")
338call network%display()
339datetime=dbadatetime(datetime_new(2014,01,06,18,00))
340call datetime%display()
341
342! can set metadata step by step
343call session%set(level=level)
344call session%set(timerange=timerange)
345call session%set(ana=ana)
346call session%set(network=network)
347call session%set(datetime=datetime)
348
349! I can use the reverse vision step by step
350!call level%dbaset(session)
351!call timerange%dbaset(session)
352!call ana%dbaset(session)
353!call network%dbaset(session)
354!call datetime%dbaset(session)
355
356! create an omogeneous vector of data
357allocate (data(2),source=[dbadatar("B12102",265.33),dbadatar("B12101",273.15)])
358
359!set and display omogeneous data
360do i =1,size(data)
361 call data(i)%display()
362 call session%set(data=data(i))
363 !call data(i)%dbaset(session)
364end do
365
366!write it in dsn or file
367call session%prendilo()
368!!$session%enq(metadata)
369
370!close message if I am writing on file
371call session%close_message()
372
373
374end subroutine write2
375
376
377subroutine write3()
378
379type(dbametadata) :: metadata
380type(dbadcv) :: datav
381type(dbalevel) :: level
382type(dbatimerange) :: timerange
383type(dbaana) :: ana
384type(dbanetwork) :: network
385type(dbadatetime) :: datetime
386
387print *,"----------------------------------------------"
388print *,"--------------- write3 ------------------------"
389
390
391!clear the dballe session
392call session%unsetall()
393
394! set metadata
395level=dbalevel(level1=105, l1=2000)
396timerange=dbatimerange(timerange=4, p1=3600,p2=7200)
397ana=dbaana(lon=12.d0,lat=45.d0)
398network=dbanetwork("generic")
399datetime=dbadatetime(datetime_new(2014,01,06,18,00))
400
401!assemble metadata
402metadata=dbametadata(level=level,timerange=timerange,ana=ana,network=network,datetime=datetime)
403call metadata%display()
404
405! I can set metadata one shot
406call session%set(metadata)
407! or in the reverse vision
408!call metadata%dbaset(session)
409
410call metadata%display()
411
412! create and display an etherogeneous ensamble of data
413allocate (datav%dcv(2))
414allocate (datav%dcv(1)%dat,source=dbadatai("B12102",26312))
415allocate (datav%dcv(2)%dat,source=dbadatar("B12101",273.15))
416call datav%display()
417!set data
418call session%set(datav=datav)
419
420! or in the reverse vision
421!call datav%dbaset(session)
422
423!write it in dsn
424call session%prendilo()
425!!$session%enq(metadata)
426
427!close message if I am writing on file
428call session%close_message()
429
430end subroutine write3
431
432
433subroutine write4()
434
435type(dbametadata) :: metadata
436type(dbadataattrv) :: dataattrv
437type(dbadcv) :: attrv
438type(dbadcv) :: datav
439
440print *,"----------------------------------------------"
441print *,"--------------- write4 ------------------------"
442
443! clear the dballe session
444call session%unsetall()
445
446! define metadata
447metadata=dbametadata( &
448 level=dbalevel(level1=105, l1=2000) &
449 ,timerange=dbatimerange(timerange=4, p1=3600,p2=7200) &
450 ,ana=dbaana(lon=13.d0,lat=45.d0) &
451 ,network=dbanetwork("generic") &
452 ,datetime=dbadatetime(datetime_new(2014,01,06,18,00)))
453
454call session%set(metadata)
455
456! create and display an etherogeneous ensamble of data
457allocate (datav%dcv(2))
458allocate (datav%dcv(1)%dat,source=dbadatai("B12102",26312))
459allocate (datav%dcv(2)%dat,source=dbadatar("B12101",273.15))
460
461! create and display an etherogeneous ensamble of attr
462allocate (attrv%dcv(3))
463allocate (attrv%dcv(1)%dat,source=dbadatai("*B33192",30))
464allocate (attrv%dcv(2)%dat,source=dbadatac("*B33193","70"))
465allocate (attrv%dcv(3)%dat,source=dbadatad("*B33194",50.d0))
466call attrv%display()
467
468! assemble data and attribute
469allocate (dataattrv%dataattr(2))
470! first with attribute
471allocate (dataattrv%dataattr(1)%dat,source=datav%dcv(1)%dat)
472dataattrv%dataattr(1)%attrv=attrv
473! second without attribute
474allocate (dataattrv%dataattr(2)%dat,source=datav%dcv(2)%dat)
475allocate (dataattrv%dataattr(2)%attrv%dcv(0))
476
477call dataattrv%display()
478
479! write data and attribute
480!call session%extrude(dataattrv)
481call dataattrv%extrude(session)
482
483
484! work on constant station data
485call session%set(metadata%dbacontextana())
486
487! write the same data and attribute as constant station data
488!call session%extrude(dataattrv)
489call dataattrv%extrude(session)
490
491end subroutine write4
492
493
494subroutine write5()
495
496print *,"----------------------------------------------"
497print *,"--------------- write5 ------------------------"
498
499call session%messages_open_output(filename="example_dballe_write5.bufr")
500
501call session%set(filter=dbafilter(coordmin=dbacoord(lon= 9.d0,lat=44.d0),&
502 coordmax=dbacoord(lon=11.d0,lat=46.d0)))
503call session%messages_write_next()
504
505end subroutine write5
506
507
508
509subroutine read1()
510
511type(dbametaanddata),allocatable :: metaanddatav(:)
512type(dbaana),allocatable :: ana(:)
513type(dbafilter) :: filter
514integer :: i
515
516print *,"----------------------------------------------"
517print *,"--------------- read1 ------------------------"
518
519call session%set(filter=dbafilter())
520
521print *, "anagrafica:"
522filter=dbafilter(coordmin=dbacoord(lon= 9.d0,lat=44.d0),&
523 coordmax=dbacoord(lon=11.d0,lat=46.d0))
524call session%set(filter=filter)
525call session%ingesta(ana)
526do i=1,size(ana)
527 call ana(i)%display()
528end do
529deallocate(ana)
530
531print *, "dati di anagrafica:"
532call session%set(filter=dbafilter(contextana=.true.))
533call session%ingest(metaanddatav)
534do i=1,size(metaanddatav)
535 call metaanddatav(i)%display()
536end do
537deallocate(metaanddatav)
538
539print *, "dati dati:"
540filter=dbafilter(var="B12102")
541call filter%display()
542call session%set(filter=filter)
543call session%ingest(metaanddatav)
544
545do i=1,size(metaanddatav)
546 call metaanddatav(i)%display()
547end do
548
549deallocate(metaanddatav)
550
551end subroutine read1
552
553subroutine read1lf()
554
555type(dbametaanddatalist) :: metaanddatal
556type(dbaanalist) :: anal
557logical :: status
558
559print *,"----------------------------------------------"
560print *,"--------------- read1lf ------------------------"
561
562! here (reading from file) we cannot use filters !!!
563
564! rewind
565call session%filerewind()
566
567print *, "anagrafica:"
568call session%ingesta(anal)
569call anal%display()
570status=anal%delete()
571
572! rewind
573call session%filerewind()
574
575print *, "dati di anagrafica:"
576call session%ingest(metaanddatal,filter=dbafilter(contextana=.true.))
577call metaanddatal%display()
578status = metaanddatal%delete()
579
580! rewind
581call session%filerewind()
582
583print *, "dati dati:"
584metaanddatal=dbametaanddatalist()
585call session%ingest(metaanddatal)
586call metaanddatal%display()
587status = metaanddatal%delete()
588
589end subroutine read1lf
590
591subroutine read2()
592
593type(dbametaanddata),allocatable :: metaanddatav(:)
594type(dbafilter) :: filter
595integer :: i
596type(dbadcv) :: vars,starvars
597integer :: ivalue
598real :: rvalue
599character(len=128) :: cvalue
600
601print *,"----------------------------------------------"
602print *,"--------------- read2 ------------------------"
603
604allocate (vars%dcv(2))
605allocate (vars%dcv(1)%dat,source=dbadatai("B12101"))
606allocate (vars%dcv(2)%dat,source=dbadatac("B12102"))
607
608filter=dbafilter(vars=vars)
609print *, "filter:"
610call filter%display()
611call session%set(filter=filter)
612call session%ingest(metaanddatav,filter=filter,noattr=.true.)
613
614print *, "dati dati:"
615do i=1,size(metaanddatav)
616 call metaanddatav(i)%display()
617end do
618
619
620!!!!!!!!!!!!!!!!
621! how use values
622associate(dato => metaanddatav(1)%dataattrv%dataattr(1)%dat)
623print *,dato%btable
624! with cast
625select type (dato)
626type is (dbadatai)
627 print *,"cast integer value",dato%value
628type is (dbadatar)
629 print *,"cast real value",dato%value
630type is (dbadatac)
631 print *,"cast character value",dato%value
632end select
633
634! calling %get
635call dato%get(ivalue)
636if (c_e(ivalue)) print *,"get integer value",ivalue
637call dato%get(rvalue)
638if (c_e(rvalue)) print *,"get real value",rvalue
639call dato%get(cvalue)
640if (c_e(cvalue)) print *,"get char value",cvalue
641end associate
642!!!!!!!!!!!!!!!!
643
644
645deallocate (vars%dcv)
646allocate (vars%dcv(2))
647allocate (vars%dcv(1)%dat,source=dbadatai("B13003"))
648allocate (vars%dcv(2)%dat,source=dbadatac("B12102"))
649allocate (starvars%dcv(3))
650allocate (starvars%dcv(1)%dat,source=dbadatai("*B33192"))
651allocate (starvars%dcv(2)%dat,source=dbadatac("*B33193"))
652allocate (starvars%dcv(3)%dat,source=dbadatad("*B33194"))
653
654filter=dbafilter(vars=vars,starvars=starvars)
655print *, "filter:"
656call filter%display()
657call session%set(filter=filter)
658call session%ingest(metaanddatav,filter=filter)
659
660print *, "dati dati:"
661do i=1,size(metaanddatav)
662 call metaanddatav(i)%display()
663end do
664
665deallocate(metaanddatav)
666
667end subroutine read2
668
669
670
671subroutine read2lf()
672
673type(dbametaanddatalist) :: metaanddatal
674type(dbafilter) :: filter
675type(dbadcv) :: vars,starvars
676logical :: status
677
678print *,"----------------------------------------------"
679print *,"--------------- read2lf ------------------------"
680
681! rewind
682call session%filerewind()
683
684allocate (vars%dcv(2))
685allocate (vars%dcv(1)%dat,source=dbadatai("B13003"))
686allocate (vars%dcv(2)%dat,source=dbadatac("B12102"))
687allocate (starvars%dcv(3))
688allocate (starvars%dcv(1)%dat,source=dbadatai("*B33192"))
689allocate (starvars%dcv(2)%dat,source=dbadatac("*B33193"))
690allocate (starvars%dcv(3)%dat,source=dbadatad("*B33194"))
691
692filter=dbafilter(vars=vars,starvars=starvars)
693!print *, "filter:"
694!call filter%display()
695
696call session%ingest(metaanddatal,filter=filter)
697
698print *, "dati dati:"
699call metaanddatal%display()
700
701status=metaanddatal%delete()
702
703end subroutine read2lf
704
705
706subroutine read3()
707
708type(dbametaanddatar),allocatable :: metaanddatarv(:)
709type(dbafilter) :: filter
710integer :: i
711
712print *,"----------------------------------------------"
713print *,"--------------- read3 ------------------------"
714
715print *, "dati dati:"
716filter=dbafilter(var="B12102")
717call filter%display()
718call session%set(filter=filter)
719call session%ingest(metaanddatarv)
720
721do i=1,size(metaanddatarv)
722 call metaanddatarv(i)%display()
723end do
724
725print *,"max=",maxval(metaanddatarv(:)%dbadatar%value)
726print *,"livelli", metaanddatarv(:)%metadata%level
727
728deallocate(metaanddatarv)
729
730end subroutine read3
731
732
733
734!!$subroutine read3l()
735!!$
736!!$type(dbametaanddatarlist) :: metaanddatarl
737!!$type(dbametaanddatar),allocatable :: metaanddatarv(:)
738!!$type(dbafilter) :: filter
739!!$logical :: status
740!!$
741!!$print *,"----------------------------------------------"
742!!$print *,"--------------- read3l ------------------------"
743!!$
744!!$! rewind
745!!$call session%filerewind()
746!!$
747!!$print *, "dati dati:"
748!!$filter=dbafilter(var="B12102")
749!!$call session%set(filter=filter)
750!!$call session%ingest(metaanddatarl)
751!!$
752!!$call metaanddatarl%display()
753!!$
754!!$metaanddatarv=metaanddatarl%toarray()
755!!$print *,"max=",maxval(metaanddatarv(:)%dbadatar%value)
756!!$
757!!$status = metaanddatarl%delete()
758!!$
759!!$end subroutine read3l
760
761
762subroutine read3lfmem()
763
764type(dbametaanddatarlist) :: metaanddatarl
765type(dbametaanddatar),allocatable :: metaanddatarv(:)
766type(dbafilter) :: filter
767logical :: status
768integer :: i
769
770print *,"----------------------------------------------"
771print *,"--------------- read3lfmem ------------------------"
772
773! connect to dsn type DBA
774connection=dbaconnection(dsn="mem:")
775session=dbasession(connection,wipe=.true.,anaflag="write", dataflag="write", attrflag="write",memdb=.true.)
776call session%messages_open_input(filename="example_dballe.bufr",mode="r",format="BUFR",simplified=.true.)
777
778filter=dbafilter(var="B12101")
779
780! data
781do while (session%messages_read_next())
782
783 call session%set(filter=filter)
784
785 call session%ingest()
786 call session%ingest(metaanddatarv)
787
788 do i=1,size(metaanddatarv)
789 call metaanddatarv(i)%display()
790 call metaanddatarl%append(metaanddatarv(i))
791 end do
792
793 call session%remove_all()
794 deallocate (metaanddatarv)
795end do
796
797
798print *, "dati dati:"
799call metaanddatarl%display()
800
801metaanddatarv=metaanddatarl%toarray()
802print *,"max=",maxval(metaanddatarv(:)%dbadatar%value)
803
804status = metaanddatarl%delete()
805
806!close everythings
807call session%delete()
808call connection%delete()
809
810end subroutine read3lfmem
811
812
813subroutine readmem()
814
815type(dbametaanddatalist) :: metaanddatal
816type(dbafilter) :: filter
817logical :: status
818
819print *,"----------------------------------------------"
820print *,"--------------- readmem ------------------------"
821
822! create mem db where put bufr from file
823session=dbasession(wipe=.true.,&
824 filename="example_dballe.bufr",mode="r",format="BUFR",memdb=.true.,loadfile=.false.)
825
826filter=dbafilter(var="B12101")
827! in this case we have to pass filter and do not set it
828call session%ingest(metaanddatal,filter=filter)
829
830print *, "dati dati:"
831call metaanddatal%display()
832
833status = metaanddatal%delete()
834
835call session%delete()
836
837end subroutine readmem
838
839
840subroutine read4l()
841
842type(dbametaanddata) :: myelement
843type(dbametaanddatalist) :: metaanddatal
844type(dbafilter) :: filter
845
846print *,"----------------------------------------------"
847print *,"--------------- read4l ------------------------"
848
849print *, "dati dati:"
850filter=dbafilter(var="B12102")
851!call session%set(filter=filter)
852
853print*,"prima di ingest"
854call session%ingest(metaanddatal,filter=filter)
855print*,"prima di display"
856call metaanddatal%display()
857
858print *,"number of list elements=",metaanddatal%countelements()
859print *,"seek return status =", metaanddatal%seek(1)
860myelement=metaanddatal%current()
861print *,"list index 1 ="
862call myelement%display()
863print *,"status delete=", metaanddatal%delete()
864
865end subroutine read4l
866
867
868subroutine delete1()
869
870type(dbafilter) :: filter
871
872print *,"----------------------------------------------"
873print *,"--------------- delete1 ----------------------"
874
875filter=dbafilter(var="B12101")
876call session%set(filter=filter)
877call session%dissolve()
878
879end subroutine delete1
880
881
882subroutine delete2()
883
884type(dbametadata),allocatable :: metadata(:)
885type(dbafilter) :: filter
886
887print *,"----------------------------------------------"
888print *,"--------------- delete2 ----------------------"
889
890
891allocate(metadata(1))
892
893metadata(1)=dbametadata( &
894 level=dbalevel(level1=105, l1=2000) &
895 ,timerange=dbatimerange(timerange=4, p1=3600,p2=7200) &
896 ,ana=dbaana(lon=11.d0,lat=45.d0) &
897 ,network=dbanetwork("generic") &
898 ,datetime=dbadatetime(datetime_new(2014,01,06,18,00)))
899
900
901filter=dbafilter(var="B12102")
902call session%set(filter=filter)
903call session%dissolve(metadata)
904
905deallocate(metadata)
906
907end subroutine delete2
908
909
910
911subroutine delete3()
912
913type(dbametadata),allocatable :: metadata(:)
914type(dbafilter) :: filter
915
916print *,"----------------------------------------------"
917print *,"--------------- delete3 ----------------------"
918
919
920allocate(metadata(1))
921
922metadata(1)=dbametadata( &
923 level=dbalevel(level1=105, l1=2000) &
924 ,timerange=dbatimerange(timerange=4, p1=3600,p2=7200) &
925 ,ana=dbaana(lon=13.d0,lat=45.d0) &
926 ,network=dbanetwork("generic") &
927 ,datetime=dbadatetime(datetime_new(2014,01,06,18,00)))
928
929filter=dbafilter(var="B12102",starvarlist="*B33194,*B33193")
930call session%set(filter=filter)
931call session%dissolveattr(metadata)
932
933deallocate(metadata)
934
935end subroutine delete3
936
937subroutine copy1()
938type(dbasession) :: sessioncp
939type(dbaconnection) :: connectioncp
940type(dbametaanddata):: metaanddata
941
942print *,"----------------------------------------------"
943print *,"--------------- copy1 ----------------------"
944
945! connect to dsn
946connectioncp=dbaconnection(dsn="sqlite:/tmp/dballecopy1.sqlite")
947!connectioncp=dbaconnection(dsn="mem:/tmp/dballecopy1")
948sessioncp=dbasession(connectioncp,wipe=.true.,anaflag="write", dataflag="write", attrflag="write")
949
950! data
951call session%set(filter=dbafilter())
952call session%ingest()
953do while (c_e(session%count) .and. session%count > 0 )
954 call session%ingest (metaanddata)
955 !call sessioncp%extrude(metaanddata)
956 call metaanddata%extrude(sessioncp)
957 if (session%file) call session%close_message()
958end do
959
960! constant data
961call session%set(filter=dbafilter(contextana=.true.))
962call session%ingest()
963do while (c_e(session%count) .and. session%count > 0)
964 call session%ingest (metaanddata)
965 !call sessioncp%extrude(metaanddata)
966 call metaanddata%extrude(sessioncp)
967 if (session%file) call session%close_message()
968end do
969
970
971!close everythings
972call sessioncp%delete()
973call connectioncp%delete()
974
975end subroutine copy1
976
977
978subroutine copy1f()
979type(dbametaanddata), allocatable:: metaanddatav(:)
980integer :: i
981
982print *,"----------------------------------------------"
983print *,"--------------- copy1f ----------------------"
984
985! data
986call sessionfrom%filerewind()
987
988call sessionfrom%set(filter=dbafilter())
989call sessionfrom%ingest(metaanddatav)
990do while (size(metaanddatav) >0)
991 print*,"read/write data; count: ",sessionfrom%count
992 do i =1,size(metaanddatav)
993 print *, "display metaanddatav index: ", i
994! call metaanddatav(i)%display()
995 !call sessionto%extrude(metaanddatav(i))
996 call metaanddatav(i)%extrude(sessionto)
997 end do
998! call sessionto%close_message()
999 call sessionfrom%ingest(metaanddatav)
1000end do
1001deallocate (metaanddatav)
1002
1003! constant data
1004call sessionfrom%filerewind()
1005
1006call sessionfrom%set(filter=dbafilter(contextana=.true.))
1007call sessionfrom%ingest(metaanddatav)
1008do while (size(metaanddatav) >0)
1009 print*,"read/write data; count: ",sessionfrom%count
1010 do i =1,size(metaanddatav)
1011 print *, "display metaanddatav index: ", i
1012! call metaanddatav(i)%display()
1013 !call sessionto%extrude(metaanddatav(i))
1014 call metaanddatav(i)%extrude(sessionto)
1015 end do
1016! call sessionto%close_message()
1017 call sessionfrom%ingest(metaanddatav)
1018end do
1019
1020end subroutine copy1f
1021
1022
1023subroutine copy1fmem()
1024type(dbametaanddata),allocatable :: metaanddatav(:)
1025integer :: i
1026
1027print *,"----------------------------------------------"
1028print *,"--------------- copy1fmem ----------------------"
1029
1030
1031do while (sessionfrom%messages_read_next())
1032 print*,"read/write message"
1033 call sessionfrom%set(filter=dbafilter())
1034 call sessionfrom%ingest(metaanddatav)
1035 do i =1,size(metaanddatav)
1036 call metaanddatav(i)%display()
1037 !call sessionto%extrude(metaanddatav(i))
1038 call metaanddatav(i)%extrude(sessionto)
1039 end do
1040 call sessionto%prendilo()
1041
1042 print *,"contextana"
1043 call sessionfrom%set(filter=dbafilter(contextana=.true.))
1044 call sessionfrom%ingest(metaanddatav)
1045 do i =1,size(metaanddatav)
1046 call metaanddatav(i)%display()
1047 !call sessionto%extrude(metaanddatav(i))
1048 call metaanddatav(i)%extrude(sessionto)
1049 end do
1050 call sessionto%close_message()
1051 call sessionfrom%remove_all()
1052end do
1053
1054end subroutine copy1fmem
1055
1056subroutine copy1f2mem()
1057
1058print *,"----------------------------------------------"
1059print *,"--------------- copy1f2mem ----------------------"
1060
1061do while (sessionfrom%messages_read_next())
1062 call sessionfrom%messages_write_next()
1063 call sessionfrom%remove_all()
1064end do
1065
1066end subroutine copy1f2mem
1067
1068
1069subroutine copy2()
1070type(dbasession) :: sessioncp
1071type(dbaconnection) :: connectioncp
1072type(dbametaanddatac):: metaanddatac
1073
1074print *,"----------------------------------------------"
1075print *,"--------------- copy2 ----------------------"
1076
1077! connect to dsn
1078connectioncp=dbaconnection(dsn="sqlite:/tmp/dballecopy2.sqlite")
1079sessioncp=dbasession(connectioncp,wipe=.true.,anaflag="write", dataflag="write", attrflag="write")
1080
1081! use character type
1082
1083! data
1084call session%set(filter=dbafilter())
1085call session%ingest_metaanddatac()
1086do while (c_e(session%count) .and. session%count > 0)
1087 call session%ingest_metaanddatac(metaanddatac)
1088 !call sessioncp%extrude(metaanddatac)
1089 call metaanddatac%extrude(sessioncp)
1090end do
1091
1092! constant data
1093call session%set(filter=dbafilter(contextana=.true.))
1094call session%ingest_metaanddatac()
1095do while (c_e(session%count) .and. session%count > 0)
1096 call session%ingest_metaanddatac(metaanddatac)
1097 !call sessioncp%extrude(metaanddatac)
1098 call metaanddatac%extrude(sessioncp)
1099end do
1100
1101
1102!close everythings
1103call sessioncp%delete()
1104call connectioncp%delete()
1105
1106end subroutine copy2
1107
1108end program example_dballe
Emit log message for a category with specific priority.
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
fortran 2003 interface to geo_coord
character version for dbadata
doubleprecision version for dbadata
integer version for dbadata
real version for dbadata
filter to apply before ingest data
double linked list of dbametaanddata
summ of all metadata pieces
manage session handle

Generated with Doxygen.