Data ListingΒΆ
Shown below is a simple program that decodes all the data in an ODB-2 file and outputs it to stdout.
/**
* To build this program, please make sure to reference linked library:
*
* gcc -lodccore -o odc-c-ls odc_ls.c
*/
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "odc/api/odc.h"
#define CHECK_RESULT(x) \
do { \
int rc = (x); \
if (rc != ODC_SUCCESS) { \
fprintf(stderr, "Error calling odc function \"%s\": %s\n", #x, odc_error_string(rc)); \
exit(1); \
} \
} while (false); \
void usage() {
fprintf(stderr, "Usage:\n");
fprintf(stderr, " odc-c-ls <odb2 file>\n\n");
}
void write_header(odc_frame_t* frame, int ncols) {
int col;
for (col = 0; col < ncols; ++col) {
const char* name;
int type;
int element_size;
int bitfield_count;
odc_frame_column_attributes(frame, col, &name, &type, &element_size, &bitfield_count);
int padding = element_size - floor(log10(abs(col + 1))) + 2;
fprintf(stdout, "%d. %-*s\t", col + 1, padding, name);
}
fprintf(stdout, "\n");
}
void write_bitfield(int64_t val, int nbits) {
int bit;
for (bit = nbits-1; bit >= 0; --bit) {
fprintf(stdout, "%s", (val & (1 << bit)) ? "1" : "0");
}
}
void write_data(odc_decoder_t* decoder, odc_frame_t* frame, long nrows, int ncols) {
const void* data;
long width;
long height;
bool columnMajor;
CHECK_RESULT(odc_decoder_data_array(decoder, &data, &width, &height, &columnMajor));
const char* name;
int bitfield_count;
int column_types[ncols];
int column_offsets[ncols];
int column_sizes[ncols];
int column_bf_sizes[ncols];
long missing_integer;
double missing_double;
int current_offset = 0;
int col;
for (col = 0; col < ncols; ++col) {
CHECK_RESULT(odc_frame_column_attributes(frame, col, &name, &column_types[col],
&column_sizes[col], &bitfield_count));
column_offsets[col] = current_offset;
current_offset += column_sizes[col];
if (column_types[col] == ODC_BITFIELD) {
int bf;
for (bf = 0; bf < bitfield_count; ++bf) {
const char* bf_name;
int bf_offset;
int bf_size;
CHECK_RESULT(odc_frame_bitfield_attributes(frame, col, bf, &bf_name, &bf_offset, &bf_size));
if (bf == 0) column_bf_sizes[col] = bf_size;
else column_bf_sizes[col] += bf_size;
}
}
}
CHECK_RESULT(odc_missing_integer(&missing_integer));
CHECK_RESULT(odc_missing_double(&missing_double));
int row;
for (row = 0; row < nrows; row++) {
int col;
for (col = 0; col < ncols; col++) {
const void* value_p = &((const char*)data)[column_offsets[col] + (width * row)];
switch (column_types[col]) {
case ODC_INTEGER:
if (*(const int64_t*)value_p == missing_integer) {
fprintf(stdout, ".");
}
else {
fprintf(stdout, "%-*lld", column_sizes[col], *(const int64_t*)value_p);
}
break;
case ODC_BITFIELD:
if (*(const int64_t*)value_p == 0) {
fprintf(stdout, ".");
}
else {
write_bitfield(*(const int64_t*)value_p, column_bf_sizes[col]);
}
break;
case ODC_REAL:
case ODC_DOUBLE:
if (*(const double*)value_p == missing_double) {
fprintf(stdout, ".");
}
else {
fprintf(stdout, "%-*f", column_sizes[col], *(const double*)value_p);
}
break;
case ODC_STRING:
fprintf(stdout, "%-.*s", column_sizes[col], (const char*)value_p);
break;
default:
fprintf(stdout, "<unknown>");
break;
}
fprintf(stdout, "\t");
}
fprintf(stdout, "\n");
}
}
int main(int argc, char* argv[]) {
if (argc != 2) {
usage();
return 1;
}
char* path = argv[1];
CHECK_RESULT(odc_initialise_api()); // initialising api
CHECK_RESULT(odc_integer_behaviour(ODC_INTEGERS_AS_LONGS)); // change from default doubles
odc_reader_t* reader = NULL;
odc_frame_t* frame = NULL;
CHECK_RESULT(odc_open_path(&reader, path)); // opening path
CHECK_RESULT(odc_new_frame(&frame, reader)); // initialising frame
int rc;
long max_aggregated_rows = 1000000;
while ((rc = odc_next_frame_aggregated(frame, max_aggregated_rows)) == ODC_SUCCESS) {
int ncols;
CHECK_RESULT(odc_frame_column_count(frame, &ncols)); // getting column count
write_header(frame, ncols);
long nrows;
odc_decoder_t* decoder = NULL;
CHECK_RESULT(odc_new_decoder(&decoder)); // initialising decoder
CHECK_RESULT(odc_decoder_defaults_from_frame(decoder, frame)); // setting decoder structure
CHECK_RESULT(odc_decode(decoder, frame, &nrows)); // decoding data
write_data(decoder, frame, nrows, ncols);
CHECK_RESULT(odc_free_decoder(decoder)); // cleaning up decoder
}
if (rc != ODC_ITERATION_COMPLETE) {
fprintf(stderr, "Error: %s\n", odc_error_string(rc)); // unsuccessful end of frame iteration
return 1;
}
CHECK_RESULT(odc_free_frame(frame));
CHECK_RESULT(odc_close(reader)); // closing reader
return 0;
}
To use this sample program, invoke it from the command line with a path to an ODB-2 file:
./odc-c-ls example.odb
1. int_column 2. real_column 3. str8_column 4. str16_column 5. bitfield_column
1 1.200000 abcdefgh abcdefghijklmnop 0000001
2 3.400000 01234567 0123456789012345 0001011
3 5.600000 string another string 1101011
/**
* To build this program, please make sure to reference linked libraries:
*
* g++ -std=c++11 -leckit -lodccore -o odc-cpp-ls odc_ls.cc
*/
#include <iostream>
#include <iomanip>
#include <cstdint>
#include <vector>
#include <numeric>
#include <math.h>
#include "eckit/runtime/Main.h"
#include "odc/api/Odb.h"
void usage() {
std::cerr << "Usage:\n odc-cpp-ls <odb2 file>" << std::endl << std::endl;
}
void write_header(int i, odc::api::ColumnInfo col) {
int padding = col.decodedSize - floor(log10(abs(i + 1))) + 2;
std::cout << (i + 1) << ". " << std::left << std::setw(padding) << col.name << "\t";
}
void write_data(size_t nrows, size_t ncols, std::vector<odc::api::ColumnInfo> columnInfo,
std::vector<odc::api::StridedData> strides, std::vector<int> nbits) {
long integer_missing = odc::api::Settings::integerMissingValue();
double double_missing = odc::api::Settings::doubleMissingValue();
for (size_t row = 0; row < nrows; ++row) {
for (size_t col = 0; col < ncols; ++col) {
switch (columnInfo[col].type) {
case odc::api::INTEGER: {
int64_t val = *reinterpret_cast<const int64_t*>(strides[col][row]);
if (val == integer_missing) {
std::cout << ".";
}
else {
std::cout << std::left << std::setw(columnInfo[col].decodedSize) << val;
}
break;
}
case odc::api::BITFIELD: {
int64_t val = *reinterpret_cast<const int64_t*>(strides[col][row]);
for (int bit = nbits[col]-1; bit >= 0; --bit) {
std::cout << ((val & (1 << bit)) ? "1" : "0");
}
break;
}
case odc::api::REAL:
case odc::api::DOUBLE: {
double val = *reinterpret_cast<const double*>(strides[col][row]);
if (val == double_missing) {
std::cout << ".";
} else {
std::cout << std::left << std::setw(columnInfo[col].decodedSize) << val;
}
break;
}
std::cout << *reinterpret_cast<const int64_t*>(strides[col][row]);
break;
case odc::api::STRING:
std::cout << std::left << std::setw(columnInfo[col].decodedSize) <<
std::string(strides[col][row],
::strnlen(strides[col][row], columnInfo[col].decodedSize));
break;
default:
ASSERT(false);
};
std::cout << "\t";
}
std::cout << std::endl;
}
}
int main(int argc, char** argv) {
if (argc != 2) {
usage();
return 1;
}
char* path = argv[1];
eckit::Main::initialise(argc, argv);
// Initialise library
odc::api::Settings::treatIntegersAsDoubles(false);
// Open the specified ODB file
odc::api::Reader reader(path);
// Iterate through the frames
odc::api::Frame frame;
while ((frame = reader.next())) {
// Properties of this frame
size_t ncols = frame.columnCount();
size_t nrows = frame.rowCount();
const auto& columnInfo = frame.columnInfo();
// Print headers & determine storage requirements
size_t row_size = 0;
size_t i;
for (i = 0; i < frame.columnCount(); ++i) {
const auto& col(columnInfo[i]);
write_header(i, col);
row_size += col.decodedSize;
}
std::cout << std::endl;
// Allocate storage required
size_t storage_size = row_size * nrows;
std::vector<char> buffer(storage_size);
// Decoder prerequisites
std::vector<std::string> columns;
std::vector<odc::api::StridedData> strides;
std::vector<int> nbits;
char* ptr = &buffer[0];
for (const auto& col : columnInfo) {
columns.push_back(col.name);
strides.emplace_back(odc::api::StridedData{ptr, nrows, col.decodedSize, col.decodedSize});
ptr += nrows * col.decodedSize;
// Bits for formatting bitfields
if (col.type == odc::api::BITFIELD) {
nbits.push_back(
std::accumulate(col.bitfield.begin(), col.bitfield.end(), 0,
[](int l, const odc::api::ColumnInfo::Bit& r) {
return l + r.size;
}));
}
else {
nbits.push_back(0);
}
}
ASSERT(ptr == (&buffer[0] + storage_size));
// Decode the data
int nthreads = 4;
odc::api::Decoder decoder(columns, strides);
decoder.decode(frame, nthreads);
// Iterate through the decoded data, and print it
write_data(nrows, ncols, columnInfo, strides, nbits);
}
return 0;
}
To use this sample program, invoke it from the command line with a path to an ODB-2 file:
./odc-cpp-ls example.odb
1. int_column 2. real_column 3. str8_column 4. str16_column 5. bitfield_column
1 1.2 abcdefgh abcdefghijklmnop 0000001
2 3.4 01234567 0123456789012345 0001011
3 5.6 string another string 1101011
! To build this program, please make sure to first compile odc Fortran module,
! and then reference linked libraries:
!
! gfortran -c ../../src/odc/api/odc.f90
! gfortran -lodccore -lfodc -o odc-fortran-ls odc_ls.f90
program odc_ls
use, intrinsic :: iso_fortran_env, only: error_unit,output_unit
use odc
implicit none
character(255) :: path
type(odc_reader) :: reader
type(odc_frame) :: frame
type(odc_decoder) :: decoder
integer(8) :: nrows
integer :: err, ncols
if (command_argument_count() /= 1) then
call usage()
stop 1
end if
call get_command_argument(1, path)
call check_call(odc_initialise_api(), "initialising api")
call check_call(reader%open_path(trim(path)), "opening path")
call check_call(frame%initialise(reader), "initialising frame")
err = frame%next()
do while (err == ODC_SUCCESS)
call check_call(frame%column_count(ncols), "getting column count")
call write_header(frame, ncols)
call check_call(decoder%initialise(), "initialising decoder")
call check_call(decoder%defaults_from_frame(frame), "setting decoder structure")
call check_call(decoder%decode(frame, nrows), "decoding data")
call write_data(decoder, frame, nrows, ncols)
call check_call(decoder%free(), "cleaning up decoder")
err = frame%next()
end do
if (err /= ODC_ITERATION_COMPLETE) call check_call(err, "get next frame")
call check_call(reader%close(), "closing reader")
contains
subroutine check_call(err, desc)
integer, intent(in) :: err
character(*), intent(in) :: desc
if (err /= ODC_SUCCESS) then
write(error_unit, *) '**** An error occurred in ODC library'
write(error_unit, *) 'Description: ', desc
write(error_unit, *) 'Error: ', odc_error_string(err)
stop 1
end if
end subroutine
subroutine usage()
write(output_unit, *) 'Usage:'
write(output_unit, *) ' odc-fortran-ls <odb2 file>'
end subroutine
subroutine write_header(frame, ncols)
type(odc_frame), intent(in) :: frame
integer, intent(in) :: ncols
character(:), allocatable :: name_str
integer :: col
do col = 1, ncols
call write_integer(output_unit, col)
call check_call(frame%column_attributes(col, name=name_str), "getting column name")
write(output_unit, '(3a)', advance='no') '. ', name_str, char(9)
end do
write(output_unit,*)
end subroutine
subroutine write_data(decoder, frame, nrows, ncols)
type(odc_decoder), intent(inout) :: decoder
type(odc_frame), intent(in) :: frame
integer(8), intent(in) :: nrows
integer, intent(in) :: ncols
integer(8) :: row
real(8), pointer :: array_data(:,:)
integer :: ncols_decoder, ncols_frame, col, current_index, bitfield_count
integer(8) :: missing_integer
real(8) :: missing_double
integer, dimension(ncols) :: types, sizes, bf_sizes, indexes
integer, target :: bf, bf_size
call check_call(decoder%data(array_data), "getting access to data")
call check_call(decoder%column_count(ncols_decoder), "getting decoder column count")
call check_call(frame%column_count(ncols_frame), "getting frame column count")
if (ncols_decoder /= ncols_frame .or. ncols_decoder /= ncols) then
write(error_unit, *) 'Something went wrong in the decode target initialisation'
stop 1
end if
current_index = 1
do col = 1, ncols
call check_call(frame%column_attributes(col, type=types(col), bitfield_count=bitfield_count), "getting column type")
call check_call(decoder%column_data_array(col, element_size_doubles=sizes(col)), "getting element size")
indexes(col) = current_index
current_index = current_index + sizes(col)
if (types(col) == ODC_BITFIELD) then
bf_sizes(col) = 0
do bf = 1, bitfield_count
call check_call(frame%bitfield_attributes(col, bf, size=bf_size), "getting bitfield size")
bf_sizes(col) = bf_sizes(col) + bf_size
end do
end if
end do
call check_call(odc_missing_integer(missing_integer), "getting missing integer")
call check_call(odc_missing_double(missing_double), "getting missing double")
do row = 1, nrows
do col = 1, ncols
select case(types(col))
case (ODC_INTEGER)
if (int(array_data(row, indexes(col))) == missing_integer) then
write(output_unit, '(a)', advance='no') '.'
else
call write_integer(output_unit, int(array_data(row, indexes(col))))
end if
case (ODC_BITFIELD)
if (int(array_data(row, indexes(col))) == 0) then
write(output_unit, '(a)', advance='no') '.'
else
call write_bitfield(output_unit, int(array_data(row, indexes(col))), bf_sizes(col))
end if
case (ODC_REAL, ODC_DOUBLE)
if (array_data(row, indexes(col)) == missing_double) then
write(output_unit, '(a)', advance='no') '.'
else
call write_double(output_unit, array_data(row, indexes(col)))
end if
case (ODC_STRING)
call write_string(output_unit, array_data(row, indexes(col):(indexes(col)+sizes(col)-1)))
case default
write(output_unit, '(a)', advance='no') '<unknown>'
end select
write(output_unit, '(a)', advance='no') char(9)
end do
write(output_unit, *)
end do
end subroutine
subroutine write_integer(iunit, i)
integer, intent(in) :: iunit, i
character(32) :: val
write(val, *) i
write(iunit, '(a)', advance='no') trim(adjustl(val))
end subroutine
subroutine write_bitfield(iunit, i, size)
integer, intent(in) :: iunit, i, size
character(32) :: padding_format, val
write(padding_format, '(a,i0,a)') '(b32.', size, ')'
write(val, padding_format) i
write(iunit, '(a)', advance='no') trim(adjustl(val))
end subroutine
subroutine write_double(iunit, r)
integer, intent(in) :: iunit
real(8), intent(in) :: r
character(32) :: val
write(val, '(f8.4)') r
write(iunit, '(a)', advance='no') trim(adjustl(val))
end subroutine
subroutine write_string(iunit, double_string)
integer, intent(in) :: iunit
real(8), intent(in), dimension(:) :: double_string
character(8*size(double_string)) :: strtmp
if (all(transfer(double_string, 1_8, size(double_string)) == 0)) then
write(iunit, '(a)', advance='no') '.'
else
write(iunit, '(a)', advance='no') trim(adjustl(transfer(double_string, strtmp)))
end if
end subroutine
end program
To use this sample program, invoke it from the command line with a path to an ODB-2 file:
./odc-fortran-ls example.odb
1. int_column 2. real_column 3. str8_column 4. str16_column 5. bitfield_column
1 1.2000 abcdefgh abcdefghijklmnop 0000001
2 3.4000 01234567 0123456789012345 0001011
3 5.6000 string another string 1101011