Bioconductor / HDF5Array

HDF5 backend for DelayedArray objects
https://bioconductor.org/packages/HDF5Array
9 stars 13 forks source link

write_block fails for SparseArraySeed #30

Closed LTLA closed 4 years ago

LTLA commented 4 years ago

Looks like HDF5Array's write_block hasn't been updated to handle SparseArraySeed:

library(Matrix)
y <- rsparsematrix(1000, 1000, 0.01)

library(HDF5Array)
as(y, "HDF5Array")
## Error in UseMethod("h5writeDataset") : 
##   no applicable method for 'h5writeDataset' applied to an object of class "c('SparseArraySeed', 'Array')"

I ended up just slapping the following inside write_block for the time being:

if (is(block, "SparseArraySeed")) {
    block <- as.array(block)
}

A better solution would be to use an indexed write, but I just couldn't get the following to work:

fhandle <- H5Fopen(x@filepath)
on.exit(H5Fclose(fhandle))

fspace <- H5Screate_simple(dim(x))
on.exit(H5Sclose(fspace), add=TRUE, after=FALSE)
indices <- nzindex(block)
H5Sselect_index(fspace,
    index=lapply(seq_len(ncol(indices)), function(i) indices[,i]))

dspace <- H5Screate_simple(length(nzdata(block)))
on.exit(H5Sclose(dspace), add=TRUE, after=FALSE)
H5Sselect_hyperslab(dspace)

dhandle <- H5Dopen(fhandle, x@name)
on.exit(H5Dclose(dhandle), add=TRUE, after=FALSE)
H5Dwrite(dhandle, buf=nzdata(block), h5spaceMem=dspace, h5spaceFile=fspace)
Session information ``` R version 4.0.0 Patched (2020-05-01 r78341) Platform: x86_64-pc-linux-gnu (64-bit) Running under: Ubuntu 18.04.4 LTS Matrix products: default BLAS: /home/luna/Software/R/R-4-0-branch-dev/lib/libRblas.so LAPACK: /home/luna/Software/R/R-4-0-branch-dev/lib/libRlapack.so locale: [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8 [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8 [7] LC_PAPER=en_US.UTF-8 LC_NAME=C [9] LC_ADDRESS=C LC_TELEPHONE=C [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C attached base packages: [1] parallel stats4 stats graphics grDevices utils datasets [8] methods base other attached packages: [1] HDF5Array_1.17.2 rhdf5_2.33.3 DelayedArray_0.15.4 [4] IRanges_2.23.10 S4Vectors_0.27.12 BiocGenerics_0.35.4 [7] matrixStats_0.56.0 Matrix_1.2-18 loaded via a namespace (and not attached): [1] compiler_4.0.0 tools_4.0.0 rhdf5filters_1.1.0 grid_4.0.0 [5] lattice_0.20-41 Rhdf5lib_1.11.2 ```
hpages commented 4 years ago

Fixed in HDF5Array 1.17.3 (commit ee51d6d9efd8e73d71f9d4498a41849a02a3c810). Just the simple fix for now. I think the more elaborated fix would need to use H5Sselect_elements() but rhdf5 does not expose it so that would need to be done in C.