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

Generated with Doxygen.