diff --git a/src/main/java/ru/itmo/ctlab/hict/hict_library/chunkedfile/MatrixQueries.java b/src/main/java/ru/itmo/ctlab/hict/hict_library/chunkedfile/MatrixQueries.java index c6beef5..d5f685a 100644 --- a/src/main/java/ru/itmo/ctlab/hict/hict_library/chunkedfile/MatrixQueries.java +++ b/src/main/java/ru/itmo/ctlab/hict/hict_library/chunkedfile/MatrixQueries.java @@ -24,6 +24,7 @@ package ru.itmo.ctlab.hict.hict_library.chunkedfile; +import ch.systemsx.cisd.hdf5.HDF5DataClass; import ch.systemsx.cisd.hdf5.IndexMap; import lombok.RequiredArgsConstructor; import lombok.extern.slf4j.Slf4j; @@ -72,7 +73,6 @@ public MatrixQueries.MatrixWithWeights getSubmatrix(final @NotNull ResolutionDes final var queryRows = (int) (endRowExcl - startRowIncl); final var queryCols = (int) (endColExcl - startColIncl); - final long[][] result = new long[queryRows][queryCols]; int deltaRow = (int) (startRow - startRowIncl); int deltaCol = (int) (startCol - startColIncl); @@ -96,39 +96,36 @@ public MatrixQueries.MatrixWithWeights getSubmatrix(final @NotNull ResolutionDes Objects.requireNonNull(dsBundle); final var blockMetaCache = new HashMap(); prefillBlockMetaCache(dsBundle, resolutionOrder, rowATUs, colATUs, blockMetaCache); - if (symmetricQuery) { - final var atuCount = rowATUs.size(); - var startDeltaCol = (int) (startCol - startColIncl); - for (int i = 0; i < atuCount; ++i) { - final var rowATU = rowATUs.get(i); - deltaCol = startDeltaCol; - final var rowCount = rowATU.getLength(); - for (int j = i; j < atuCount; ++j) { - final var colATU = colATUs.get(j); - final var colCount = colATU.getLength(); - fillATUIntersectionIntoWithBundle(dsBundle, resolutionDescriptor, rowATU, colATU, result, deltaRow, deltaCol, false, blockMetaCache); - if (i != j) { - for (int k = 0; k < rowCount; k++) { - for (int l = 0; l < colCount; l++) { - result[deltaCol + l][deltaRow + k] = result[deltaRow + k][deltaCol + l]; - } - } - } - deltaCol += colCount; - } - startDeltaCol += colATUs.get(i).getLength(); - deltaRow += rowCount; - } - } else { - for (final var rowATU : rowATUs) { - deltaCol = (int) (startCol - startColIncl); - for (final var colATU : colATUs) { - fillATUIntersectionIntoWithBundle(dsBundle, resolutionDescriptor, rowATU, colATU, result, deltaRow, deltaCol, false, blockMetaCache); - deltaCol += colATU.getLength(); - } - deltaRow += rowATU.getLength(); - } - } + final var result = isFloatingPointDataset(dsBundle.getReader(), dsBundle.getBlockValuesDataSet()) + ? fillSubmatrixAsDoubles( + dsBundle, + resolutionDescriptor, + rowATUs, + colATUs, + queryRows, + queryCols, + symmetricQuery, + startRow, + startRowIncl, + startCol, + startColIncl, + blockMetaCache + ) + : fillSubmatrixAsLongs( + dsBundle, + resolutionDescriptor, + rowATUs, + colATUs, + queryRows, + queryCols, + symmetricQuery, + startRow, + startRowIncl, + startCol, + startColIncl, + blockMetaCache + ); + return new MatrixQueries.MatrixWithWeights(result, paddedRowWeights, paddedColWeights, startRow, startCol, endRow, endCol, units, resolutionDescriptor); } catch (final Exception ex) { throw new RuntimeException("Dense matrix fetch failed", ex); } finally { @@ -143,7 +140,113 @@ public MatrixQueries.MatrixWithWeights getSubmatrix(final @NotNull ResolutionDes } - return new MatrixQueries.MatrixWithWeights(result, paddedRowWeights, paddedColWeights, startRow, startCol, endRow, endCol, units, resolutionDescriptor); + return new MatrixQueries.MatrixWithWeights(new LongMatrix(new long[queryRows][queryCols]), paddedRowWeights, paddedColWeights, startRow, startCol, endRow, endCol, units, resolutionDescriptor); + } + + private @NotNull RawMatrix fillSubmatrixAsLongs( + final @NotNull HDF5FileDatasetsBundle dsBundle, + final @NotNull ResolutionDescriptor resolutionDescriptor, + final @NotNull List<@NotNull ATUDescriptor> rowATUs, + final @NotNull List<@NotNull ATUDescriptor> colATUs, + final int queryRows, + final int queryCols, + final boolean symmetricQuery, + final long startRow, + final long startRowIncl, + final long startCol, + final long startColIncl, + final @NotNull Map blockMetaCache + ) { + final var result = new long[queryRows][queryCols]; + int deltaRow = (int) (startRow - startRowIncl); + int deltaCol; + if (symmetricQuery) { + final var atuCount = rowATUs.size(); + var startDeltaCol = (int) (startCol - startColIncl); + for (int i = 0; i < atuCount; ++i) { + final var rowATU = rowATUs.get(i); + deltaCol = startDeltaCol; + final var rowCount = rowATU.getLength(); + for (int j = i; j < atuCount; ++j) { + final var colATU = colATUs.get(j); + final var colCount = colATU.getLength(); + fillATUIntersectionIntoWithBundle(dsBundle, resolutionDescriptor, rowATU, colATU, result, deltaRow, deltaCol, false, blockMetaCache); + if (i != j) { + for (int k = 0; k < rowCount; k++) { + for (int l = 0; l < colCount; l++) { + result[deltaCol + l][deltaRow + k] = result[deltaRow + k][deltaCol + l]; + } + } + } + deltaCol += colCount; + } + startDeltaCol += colATUs.get(i).getLength(); + deltaRow += rowCount; + } + } else { + for (final var rowATU : rowATUs) { + deltaCol = (int) (startCol - startColIncl); + for (final var colATU : colATUs) { + fillATUIntersectionIntoWithBundle(dsBundle, resolutionDescriptor, rowATU, colATU, result, deltaRow, deltaCol, false, blockMetaCache); + deltaCol += colATU.getLength(); + } + deltaRow += rowATU.getLength(); + } + } + return new LongMatrix(result); + } + + private @NotNull RawMatrix fillSubmatrixAsDoubles( + final @NotNull HDF5FileDatasetsBundle dsBundle, + final @NotNull ResolutionDescriptor resolutionDescriptor, + final @NotNull List<@NotNull ATUDescriptor> rowATUs, + final @NotNull List<@NotNull ATUDescriptor> colATUs, + final int queryRows, + final int queryCols, + final boolean symmetricQuery, + final long startRow, + final long startRowIncl, + final long startCol, + final long startColIncl, + final @NotNull Map blockMetaCache + ) { + final var result = new double[queryRows][queryCols]; + int deltaRow = (int) (startRow - startRowIncl); + int deltaCol; + if (symmetricQuery) { + final var atuCount = rowATUs.size(); + var startDeltaCol = (int) (startCol - startColIncl); + for (int i = 0; i < atuCount; ++i) { + final var rowATU = rowATUs.get(i); + deltaCol = startDeltaCol; + final var rowCount = rowATU.getLength(); + for (int j = i; j < atuCount; ++j) { + final var colATU = colATUs.get(j); + final var colCount = colATU.getLength(); + fillATUIntersectionIntoWithBundle(dsBundle, resolutionDescriptor, rowATU, colATU, result, deltaRow, deltaCol, false, blockMetaCache); + if (i != j) { + for (int k = 0; k < rowCount; k++) { + for (int l = 0; l < colCount; l++) { + result[deltaCol + l][deltaRow + k] = result[deltaRow + k][deltaCol + l]; + } + } + } + deltaCol += colCount; + } + startDeltaCol += colATUs.get(i).getLength(); + deltaRow += rowCount; + } + } else { + for (final var rowATU : rowATUs) { + deltaCol = (int) (startCol - startColIncl); + for (final var colATU : colATUs) { + fillATUIntersectionIntoWithBundle(dsBundle, resolutionDescriptor, rowATU, colATU, result, deltaRow, deltaCol, false, blockMetaCache); + deltaCol += colATU.getLength(); + } + deltaRow += rowATU.getLength(); + } + } + return new DoubleMatrix(result); } private double @NotNull [] flattenWeights(final @NotNull List<@NotNull ATUDescriptor> atus) { @@ -402,11 +505,11 @@ public List getATUsForRange(final @NotNull ResolutionDescriptor r return reducedATUs; } - public long @NotNull [][] getATUIntersection(final @NotNull ResolutionDescriptor resolutionDescriptor, final @NotNull ATUDescriptor rowATU, final @NotNull ATUDescriptor colATU) { + public double @NotNull [][] getATUIntersection(final @NotNull ResolutionDescriptor resolutionDescriptor, final @NotNull ATUDescriptor rowATU, final @NotNull ATUDescriptor colATU) { return getATUIntersection(resolutionDescriptor, rowATU, colATU, false); } - public long @NotNull [][] getATUIntersection(final @NotNull ResolutionDescriptor resolutionDescriptor, final @NotNull ATUDescriptor rowATU, final @NotNull ATUDescriptor colATU, final boolean needsTranspose) { + public double @NotNull [][] getATUIntersection(final @NotNull ResolutionDescriptor resolutionDescriptor, final @NotNull ATUDescriptor rowATU, final @NotNull ATUDescriptor colATU, final boolean needsTranspose) { final var resolutionOrder = resolutionDescriptor.getResolutionOrderInArray(); final @NotNull var pool = this.chunkedFile.getDatasetBundlePools().get(resolutionOrder); @Nullable HDF5FileDatasetsBundle dsBundle = null; @@ -427,7 +530,7 @@ public List getATUsForRange(final @NotNull ResolutionDescriptor r } } - private long @NotNull [][] getATUIntersectionWithBundle( + private double @NotNull [][] getATUIntersectionWithBundle( final @NotNull HDF5FileDatasetsBundle dsBundle, final @NotNull ResolutionDescriptor resolutionDescriptor, final @NotNull ATUDescriptor rowATU, @@ -435,7 +538,7 @@ public List getATUsForRange(final @NotNull ResolutionDescriptor r final boolean needsTranspose, final @Nullable Map blockMetaCache ) { - final var denseMatrix = new long[needsTranspose ? colATU.getLength() : rowATU.getLength()][needsTranspose ? rowATU.getLength() : colATU.getLength()]; + final var denseMatrix = new double[needsTranspose ? colATU.getLength() : rowATU.getLength()][needsTranspose ? rowATU.getLength() : colATU.getLength()]; fillATUIntersectionIntoWithBundle(dsBundle, resolutionDescriptor, rowATU, colATU, denseMatrix, 0, 0, needsTranspose, blockMetaCache); return denseMatrix; } @@ -445,7 +548,7 @@ private void fillATUIntersectionIntoWithBundle( final @NotNull ResolutionDescriptor resolutionDescriptor, final @NotNull ATUDescriptor rowATU, final @NotNull ATUDescriptor colATU, - final long[][] target, + final double[][] target, final int dstRow, final int dstCol, final boolean needsTranspose, @@ -506,14 +609,14 @@ private void fillATUIntersectionIntoWithBundle( final var cachedPayload = this.blockDataCache.get(cacheKey); if (savedAsSparse) { - final SparseBlockData sparseBlockData; - if (cachedPayload instanceof SparseBlockData cachedSparse) { + final SparseDoubleBlockData sparseBlockData; + if (cachedPayload instanceof SparseDoubleBlockData cachedSparse) { sparseBlockData = cachedSparse; } else { - sparseBlockData = new SparseBlockData( + sparseBlockData = new SparseDoubleBlockData( reader.int64().readArrayBlockWithOffset(dsBundle.getBlockRowsDataSet(), (int) blockLength, blockOffset), reader.int64().readArrayBlockWithOffset(dsBundle.getBlockColsDataSet(), (int) blockLength, blockOffset), - reader.int64().readArrayBlockWithOffset(dsBundle.getBlockValuesDataSet(), (int) blockLength, blockOffset) + readNumericArrayBlockWithOffset(reader, dsBundle.getBlockValuesDataSet(), (int) blockLength, blockOffset) ); this.blockDataCache.put(cacheKey, sparseBlockData); } @@ -556,17 +659,12 @@ private void fillATUIntersectionIntoWithBundle( } } } else { - final DenseBlockData denseBlockData; - if (cachedPayload instanceof DenseBlockData cachedDense) { + final DenseDoubleBlockData denseBlockData; + if (cachedPayload instanceof DenseDoubleBlockData cachedDense) { denseBlockData = cachedDense; } else { final var idx = new IndexMap().bind(0, -(blockOffset + 1L)).bind(1, 0L); - final var block = reader.int64().readSlicedMDArrayBlockWithOffset( - dsBundle.getDenseBlockDataSet(), - new int[]{this.chunkedFile.getDenseBlockSize(), this.chunkedFile.getDenseBlockSize()}, - new long[]{0L, 0L}, - idx - ); + final var block = readNumericDenseBlock(reader, dsBundle.getDenseBlockDataSet(), this.chunkedFile.getDenseBlockSize(), idx); final var denseBlock = block.toMatrix(); if (blockOnMainDiagonal) { for (int i = 0; i < denseBlock.length; ++i) { @@ -575,7 +673,7 @@ private void fillATUIntersectionIntoWithBundle( } } } - denseBlockData = new DenseBlockData(denseBlock); + denseBlockData = new DenseDoubleBlockData(denseBlock); this.blockDataCache.put(cacheKey, denseBlockData); } final var denseBlock = denseBlockData.matrix(); @@ -610,6 +708,201 @@ private void fillATUIntersectionIntoWithBundle( } } + private void fillATUIntersectionIntoWithBundle( + final @NotNull HDF5FileDatasetsBundle dsBundle, + final @NotNull ResolutionDescriptor resolutionDescriptor, + final @NotNull ATUDescriptor rowATU, + final @NotNull ATUDescriptor colATU, + final long[][] target, + final int dstRow, + final int dstCol, + final boolean needsTranspose, + final @Nullable Map blockMetaCache + ) { + if (rowATU.getStripeDescriptor().stripeId() > colATU.getStripeDescriptor().stripeId()) { + fillATUIntersectionIntoWithBundle(dsBundle, resolutionDescriptor, colATU, rowATU, target, dstRow, dstCol, !needsTranspose, blockMetaCache); + return; + } + + final var resolutionOrder = resolutionDescriptor.getResolutionOrderInArray(); + final var rowStripe = rowATU.getStripeDescriptor(); + final var colStripe = colATU.getStripeDescriptor(); + final var rowStripeId = rowStripe.stripeId(); + final var colStripeId = colStripe.stripeId(); + final var blockOnMainDiagonal = (rowStripeId == colStripeId); + + final long blockIndexInDatasets = ((long) rowStripeId) * this.chunkedFile.getStripeCount()[resolutionOrder] + colStripeId; + final long blockLength; + final long blockOffset; + try { + final var reader = dsBundle.getReader(); + final BlockMeta blockMeta; + if (blockMetaCache != null) { + blockMeta = blockMetaCache.computeIfAbsent(blockIndexInDatasets, key -> { + final var blockLengthDataset = dsBundle.getBlockLengthDataSet(); + final var blockOffsetDataset = dsBundle.getBlockOffsetDataSet(); + final long[] blockLengthBuf = reader.int64().readArrayBlockWithOffset(blockLengthDataset, 1, key.longValue()); + final long[] blockOffsetBuf = reader.int64().readArrayBlockWithOffset(blockOffsetDataset, 1, key.longValue()); + return new BlockMeta(blockLengthBuf[0], blockOffsetBuf[0]); + }); + } else { + final var blockLengthDataset = dsBundle.getBlockLengthDataSet(); + final var blockOffsetDataset = dsBundle.getBlockOffsetDataSet(); + final long[] blockLengthBuf = reader.int64().readArrayBlockWithOffset(blockLengthDataset, 1, blockIndexInDatasets); + final long[] blockOffsetBuf = reader.int64().readArrayBlockWithOffset(blockOffsetDataset, 1, blockIndexInDatasets); + blockMeta = new BlockMeta(blockLengthBuf[0], blockOffsetBuf[0]); + } + blockLength = blockMeta.blockLength(); + blockOffset = blockMeta.blockOffset(); + + if (blockLength == 0L) { + return; + } + + final var firstRow = rowATU.getStartIndexInStripeIncl(); + final var firstCol = colATU.getStartIndexInStripeIncl(); + final var lastRow = rowATU.getEndIndexInStripeExcl(); + final var lastCol = colATU.getEndIndexInStripeExcl(); + final var queryRowCount = lastRow - firstRow; + final var queryColCount = lastCol - firstCol; + final var flipRows = ATUDirection.REVERSED.equals(rowATU.getDirection()); + final var flipCols = ATUDirection.REVERSED.equals(colATU.getDirection()); + final var savedAsSparse = (blockOffset >= 0L); + final var cacheKey = new BlockCacheKey(resolutionOrder, blockIndexInDatasets); + final var cachedPayload = this.blockDataCache.get(cacheKey); + + if (savedAsSparse) { + final SparseLongBlockData sparseBlockData; + if (cachedPayload instanceof SparseLongBlockData cachedSparse) { + sparseBlockData = cachedSparse; + } else { + sparseBlockData = new SparseLongBlockData( + reader.int64().readArrayBlockWithOffset(dsBundle.getBlockRowsDataSet(), (int) blockLength, blockOffset), + reader.int64().readArrayBlockWithOffset(dsBundle.getBlockColsDataSet(), (int) blockLength, blockOffset), + reader.int64().readArrayBlockWithOffset(dsBundle.getBlockValuesDataSet(), (int) blockLength, blockOffset) + ); + this.blockDataCache.put(cacheKey, sparseBlockData); + } + + for (int i = 0; i < sparseBlockData.values().length; i++) { + final var value = sparseBlockData.values()[i]; + writeSparseValueToTarget( + target, + dstRow, + dstCol, + sparseBlockData.rows()[i], + sparseBlockData.cols()[i], + value, + firstRow, + firstCol, + queryRowCount, + queryColCount, + flipRows, + flipCols, + needsTranspose + ); + if (blockOnMainDiagonal && sparseBlockData.rows()[i] != sparseBlockData.cols()[i]) { + writeSparseValueToTarget( + target, + dstRow, + dstCol, + sparseBlockData.cols()[i], + sparseBlockData.rows()[i], + value, + firstRow, + firstCol, + queryRowCount, + queryColCount, + flipRows, + flipCols, + needsTranspose + ); + } + } + } else { + final DenseLongBlockData denseBlockData; + if (cachedPayload instanceof DenseLongBlockData cachedDense) { + denseBlockData = cachedDense; + } else { + final var idx = new IndexMap().bind(0, -(blockOffset + 1L)).bind(1, 0L); + final var block = reader.int64().readSlicedMDArrayBlockWithOffset( + dsBundle.getDenseBlockDataSet(), + new int[]{this.chunkedFile.getDenseBlockSize(), this.chunkedFile.getDenseBlockSize()}, + new long[]{0L, 0L}, + idx + ); + final var denseBlock = block.toMatrix(); + if (blockOnMainDiagonal) { + for (int i = 0; i < denseBlock.length; ++i) { + for (int j = 1 + i; j < denseBlock.length; ++j) { + denseBlock[j][i] = denseBlock[i][j]; + } + } + } + denseBlockData = new DenseLongBlockData(denseBlock); + this.blockDataCache.put(cacheKey, denseBlockData); + } + final var denseBlock = denseBlockData.matrix(); + if (blockOnMainDiagonal) { + for (int outRow = 0; outRow < queryRowCount; outRow++) { + final int srcRow = flipRows ? (lastRow - 1 - outRow) : (firstRow + outRow); + for (int outCol = 0; outCol < queryColCount; outCol++) { + final int srcCol = flipCols ? (lastCol - 1 - outCol) : (firstCol + outCol); + if (needsTranspose) { + target[dstRow + outCol][dstCol + outRow] = denseBlock[srcRow][srcCol]; + } else { + target[dstRow + outRow][dstCol + outCol] = denseBlock[srcRow][srcCol]; + } + } + } + } else { + for (int outRow = 0; outRow < queryRowCount; outRow++) { + final int srcRow = flipRows ? (queryRowCount - 1 - outRow) : outRow; + for (int outCol = 0; outCol < queryColCount; outCol++) { + final int srcCol = flipCols ? (queryColCount - 1 - outCol) : outCol; + if (needsTranspose) { + target[dstRow + outCol][dstCol + outRow] = denseBlock[firstRow + srcRow][firstCol + srcCol]; + } else { + target[dstRow + outRow][dstCol + outCol] = denseBlock[firstRow + srcRow][firstCol + srcCol]; + } + } + } + } + } + } catch (Exception e) { + throw new RuntimeException(e); + } + } + + private static void writeSparseValueToTarget( + final double[][] target, + final int dstRow, + final int dstCol, + final long row, + final long col, + final double value, + final int firstRow, + final int firstCol, + final int queryRowCount, + final int queryColCount, + final boolean flipRows, + final boolean flipCols, + final boolean needsTranspose + ) { + final int localRow = (int) row - firstRow; + final int localCol = (int) col - firstCol; + if (localRow < 0 || localRow >= queryRowCount || localCol < 0 || localCol >= queryColCount) { + return; + } + final int transformedRow = flipRows ? (queryRowCount - 1 - localRow) : localRow; + final int transformedCol = flipCols ? (queryColCount - 1 - localCol) : localCol; + if (needsTranspose) { + target[dstRow + transformedCol][dstCol + transformedRow] = value; + } else { + target[dstRow + transformedRow][dstCol + transformedCol] = value; + } + } + private static void writeSparseValueToTarget( final long[][] target, final int dstRow, @@ -686,18 +979,35 @@ private record BlockMeta(long blockLength, long blockOffset) { private record BlockCacheKey(int resolutionOrder, long blockIndex) { } - private sealed interface BlockPayload permits SparseBlockData, DenseBlockData { + private sealed interface BlockPayload permits SparseDoubleBlockData, SparseLongBlockData, DenseDoubleBlockData, DenseLongBlockData { long estimatedBytes(); } - private record SparseBlockData(long[] rows, long[] cols, long[] values) implements BlockPayload { + private record SparseDoubleBlockData(long[] rows, long[] cols, double[] values) implements BlockPayload { @Override public long estimatedBytes() { - return ((long) rows.length + cols.length + values.length) * Long.BYTES; + return (((long) rows.length + cols.length) * Long.BYTES) + (((long) values.length) * Double.BYTES); } } - private record DenseBlockData(long[][] matrix) implements BlockPayload { + private record SparseLongBlockData(long[] rows, long[] cols, long[] values) implements BlockPayload { + @Override + public long estimatedBytes() { + return (((long) rows.length + cols.length + values.length) * Long.BYTES); + } + } + + private record DenseDoubleBlockData(double[][] matrix) implements BlockPayload { + @Override + public long estimatedBytes() { + if (matrix.length == 0) { + return 0L; + } + return ((long) matrix.length) * matrix[0].length * Double.BYTES; + } + } + + private record DenseLongBlockData(long[][] matrix) implements BlockPayload { @Override public long estimatedBytes() { if (matrix.length == 0) { @@ -710,6 +1020,58 @@ public long estimatedBytes() { private record BlockMetaRow(long[] blockLengths, long[] blockOffsets) { } + private static boolean isFloatingPointDataset( + final @NotNull ch.systemsx.cisd.hdf5.IHDF5Reader reader, + final @NotNull ch.systemsx.cisd.hdf5.HDF5DataSet dataSet + ) { + return reader.object().getDataSetInformation(dataSet.getDataSetPath()).getTypeInformation().getDataClass() == HDF5DataClass.FLOAT; + } + + private static double @NotNull [] readNumericArrayBlockWithOffset( + final @NotNull ch.systemsx.cisd.hdf5.IHDF5Reader reader, + final @NotNull ch.systemsx.cisd.hdf5.HDF5DataSet dataSet, + final int length, + final long offset + ) { + if (isFloatingPointDataset(reader, dataSet)) { + return reader.float64().readArrayBlockWithOffset(dataSet, length, offset); + } + final var values = reader.int64().readArrayBlockWithOffset(dataSet, length, offset); + final var result = new double[values.length]; + for (int i = 0; i < values.length; i++) { + result[i] = values[i]; + } + return result; + } + + private static @NotNull ch.systemsx.cisd.base.mdarray.MDDoubleArray readNumericDenseBlock( + final @NotNull ch.systemsx.cisd.hdf5.IHDF5Reader reader, + final @NotNull ch.systemsx.cisd.hdf5.HDF5DataSet dataSet, + final int denseBlockSize, + final @NotNull IndexMap idx + ) { + if (isFloatingPointDataset(reader, dataSet)) { + return reader.float64().readSlicedMDArrayBlockWithOffset( + dataSet, + new int[]{denseBlockSize, denseBlockSize}, + new long[]{0L, 0L}, + idx + ); + } + final var block = reader.int64().readSlicedMDArrayBlockWithOffset( + dataSet, + new int[]{denseBlockSize, denseBlockSize}, + new long[]{0L, 0L}, + idx + ); + final var longs = block.getAsFlatArray(); + final var doubles = new double[longs.length]; + for (int i = 0; i < longs.length; i++) { + doubles[i] = longs[i]; + } + return new ch.systemsx.cisd.base.mdarray.MDDoubleArray(doubles, new int[]{denseBlockSize, denseBlockSize}); + } + private static final class BlockMetaRowCache { private final int maxRows; private final LinkedHashMap cache; @@ -778,7 +1140,61 @@ private synchronized void put(final @NotNull BlockCacheKey key, final @NotNull B } } - public record MatrixWithWeights(long @NotNull [][] matrix, double @NotNull [] rowWeights, + public sealed interface RawMatrix permits LongMatrix, DoubleMatrix { + int rows(); + + int cols(); + + boolean isFloatingPoint(); + + double getAsDouble(int row, int col); + } + + public record LongMatrix(long @NotNull [][] values) implements RawMatrix { + @Override + public int rows() { + return values.length; + } + + @Override + public int cols() { + return values.length > 0 ? values[0].length : 0; + } + + @Override + public boolean isFloatingPoint() { + return false; + } + + @Override + public double getAsDouble(final int row, final int col) { + return values[row][col]; + } + } + + public record DoubleMatrix(double @NotNull [][] values) implements RawMatrix { + @Override + public int rows() { + return values.length; + } + + @Override + public int cols() { + return values.length > 0 ? values[0].length : 0; + } + + @Override + public boolean isFloatingPoint() { + return true; + } + + @Override + public double getAsDouble(final int row, final int col) { + return values[row][col]; + } + } + + public record MatrixWithWeights(@NotNull RawMatrix matrix, double @NotNull [] rowWeights, double @NotNull [] colWeights, long startRowIncl, long startColIncl, long endRowExcl, long endColExcl, @NotNull QueryLengthUnit units, @NotNull ResolutionDescriptor resolutionDescriptor) { diff --git a/src/main/java/ru/itmo/ctlab/hict/hict_library/converters/McoolToHictConverter.java b/src/main/java/ru/itmo/ctlab/hict/hict_library/converters/McoolToHictConverter.java index 26c9208..efbd277 100644 --- a/src/main/java/ru/itmo/ctlab/hict/hict_library/converters/McoolToHictConverter.java +++ b/src/main/java/ru/itmo/ctlab/hict/hict_library/converters/McoolToHictConverter.java @@ -1,7 +1,9 @@ package ru.itmo.ctlab.hict.hict_library.converters; +import ch.systemsx.cisd.base.mdarray.MDDoubleArray; import ch.systemsx.cisd.base.mdarray.MDLongArray; import ch.systemsx.cisd.hdf5.HDF5Factory; +import ch.systemsx.cisd.hdf5.HDF5DataClass; import ch.systemsx.cisd.hdf5.HDF5FloatStorageFeatures; import ch.systemsx.cisd.hdf5.HDF5IntStorageFeatures; import ch.systemsx.cisd.hdf5.IHDF5Reader; @@ -166,6 +168,8 @@ private static void writeResolutionDirect( dst.int64().setAttr(treapRoot, "stripes_count", stripeCount); final var allRowsStartIndices = src.int64().readArray(resolutionRoot + "/indexes/bin1_offset"); + final var sourceCountsPath = resolutionRoot + "/pixels/count"; + final var floatingPointSignal = isFloatingPointDataset(src, sourceCountsPath); logConsumer.accept("Resolution " + resolution + ": counting sparse and dense blocks"); final var countingProgress = new PhaseProgressTracker( @@ -174,7 +178,7 @@ private static void writeResolutionDirect( logConsumer ); final int stripeWorkers = Math.max(1, Math.min(stripeWorkersRequested, Math.max(1, stripeCount))); - final var counts = countDenseAndSparse(inputPath, resolution, stripeCount, allRowsStartIndices, stripeWorkers, countingProgress::report); + final var counts = countDenseAndSparse(inputPath, resolution, stripeCount, allRowsStartIndices, stripeWorkers, floatingPointSignal, countingProgress::report); countingProgress.finish(); final var denseBlockCount = counts.denseTotal(); logConsumer.accept("Resolution " + resolution + ": finished counting blocks, denseBlocks=" + denseBlockCount); @@ -188,18 +192,21 @@ private static void writeResolutionDirect( dst.int64().createArray(blockRowsPath, nonzeroPixelCount, safeChunkLen(nonzeroPixelCount, chunkSize), intStorageFeatures); dst.int64().createArray(blockColsPath, nonzeroPixelCount, safeChunkLen(nonzeroPixelCount, chunkSize), intStorageFeatures); - dst.int64().createArray(blockValsPath, nonzeroPixelCount, safeChunkLen(nonzeroPixelCount, chunkSize), intStorageFeatures); + createNumericArray(dst, blockValsPath, nonzeroPixelCount, chunkSize, floatingPointSignal, intStorageFeatures, floatStorageFeatures); final var totalBlockCount = (long) stripeCount * stripeCount; dst.int64().createArray(blockOffsetPath, totalBlockCount, safeChunkLen(totalBlockCount, chunkSize), intStorageFeatures); dst.int64().createArray(blockLengthPath, totalBlockCount, safeChunkLen(totalBlockCount, chunkSize), intStorageFeatures); final var denseDatasetSize = Math.max(1L, denseBlockCount); - dst.int64().createMDArray( + createNumericMDArray( + dst, denseBlocksPath, new long[]{denseDatasetSize, 1L, SUBMATRIX_SIZE, SUBMATRIX_SIZE}, new int[]{1, 1, SUBMATRIX_SIZE, SUBMATRIX_SIZE}, - intStorageFeatures + floatingPointSignal, + intStorageFeatures, + floatStorageFeatures ); final var stripeSparseBase = new long[stripeCount]; @@ -253,7 +260,8 @@ private static void writeResolutionDirect( blockValsPath, blockOffsetPath, blockLengthPath, - denseBlocksPath + denseBlocksPath, + floatingPointSignal ); } written++; @@ -288,7 +296,7 @@ private static void writeResolutionDirect( final int end = Math.min(stripeCount, batchStart + batchSize); futures.add(stripeExecutor.submit(() -> { try { - final var batchBlocks = readStripeBatch(readerHolder.get(), resolution, start, end, allRowsStartIndices); + final var batchBlocks = readStripeBatch(readerHolder.get(), resolution, start, end, allRowsStartIndices, floatingPointSignal); for (int i = 0; i < batchBlocks.length; i++) { checkInterrupted(); final int stripeIdx = start + i; @@ -354,7 +362,8 @@ private static SaveBlockResult saveStripeBlocks( final @NotNull String blockValsPath, final @NotNull String blockOffsetPath, final @NotNull String blockLengthPath, - final @NotNull String denseBlocksPath + final @NotNull String denseBlocksPath, + final boolean floatingPointSignal ) { long sparseOffset = stripeSparseOffset; long denseOffset = stripeDenseOffset; @@ -364,9 +373,6 @@ private static SaveBlockResult saveStripeBlocks( final var denseFlags = blocks.denseFlags(); final var sparseRows = blocks.sparseRows(); final var sparseCols = blocks.sparseCols(); - final var sparseVals = blocks.sparseVals(); - final var denseFlats = blocks.denseFlats(); - long stripeSparseLen = 0L; for (int i = 0; i < colStripes.length; i++) { if (!denseFlags[i]) { @@ -378,45 +384,89 @@ private static SaveBlockResult saveStripeBlocks( final long[] lengthRow = new long[stripeCount]; final long[] stripeRows = stripeSparseLen > 0 ? new long[(int) stripeSparseLen] : new long[0]; final long[] stripeCols = stripeSparseLen > 0 ? new long[(int) stripeSparseLen] : new long[0]; - final long[] stripeVals = stripeSparseLen > 0 ? new long[(int) stripeSparseLen] : new long[0]; int stripePos = 0; - - for (int i = 0; i < colStripes.length; i++) { - final int blockLen = lengths[i]; - if (blockLen <= 0) { - continue; + if (floatingPointSignal) { + final var typedBlocks = (DoubleStripeBlocks) blocks; + final var sparseVals = typedBlocks.sparseVals(); + final var denseFlats = typedBlocks.denseFlats(); + final double[] stripeVals = stripeSparseLen > 0 ? new double[(int) stripeSparseLen] : new double[0]; + for (int i = 0; i < colStripes.length; i++) { + final int blockLen = lengths[i]; + if (blockLen <= 0) { + continue; + } + final int colStripe = colStripes[i]; + + if (denseFlags[i]) { + offsetRow[colStripe] = -denseOffset - 1L; + lengthRow[colStripe] = blockLen; + final var denseMd = new MDDoubleArray(denseFlats[i], new int[]{1, 1, SUBMATRIX_SIZE, SUBMATRIX_SIZE}); + dst.float64().writeMDArrayBlockWithOffset(denseBlocksPath, denseMd, new long[]{denseOffset, 0L, 0L, 0L}); + denseOffset++; + } else { + offsetRow[colStripe] = sparseOffset; + lengthRow[colStripe] = blockLen; + final var blockRows = sparseRows[i]; + final var blockCols = sparseCols[i]; + final var blockVals = sparseVals[i]; + System.arraycopy(blockRows, 0, stripeRows, stripePos, blockLen); + System.arraycopy(blockCols, 0, stripeCols, stripePos, blockLen); + System.arraycopy(blockVals, 0, stripeVals, stripePos, blockLen); + stripePos += blockLen; + sparseOffset += blockLen; + } } - final int colStripe = colStripes[i]; - - if (denseFlags[i]) { - offsetRow[colStripe] = -denseOffset - 1L; - lengthRow[colStripe] = blockLen; - final var denseMd = new MDLongArray(denseFlats[i], new int[]{1, 1, SUBMATRIX_SIZE, SUBMATRIX_SIZE}); - dst.int64().writeMDArrayBlockWithOffset(denseBlocksPath, denseMd, new long[]{denseOffset, 0L, 0L, 0L}); - denseOffset++; - } else { - offsetRow[colStripe] = sparseOffset; - lengthRow[colStripe] = blockLen; - final var blockRows = sparseRows[i]; - final var blockCols = sparseCols[i]; - final var blockVals = sparseVals[i]; - System.arraycopy(blockRows, 0, stripeRows, stripePos, blockLen); - System.arraycopy(blockCols, 0, stripeCols, stripePos, blockLen); - System.arraycopy(blockVals, 0, stripeVals, stripePos, blockLen); - stripePos += blockLen; - sparseOffset += blockLen; + final long rowBase = (long) rowStripe * stripeCount; + dst.int64().writeArrayBlockWithOffset(blockOffsetPath, offsetRow, stripeCount, rowBase); + dst.int64().writeArrayBlockWithOffset(blockLengthPath, lengthRow, stripeCount, rowBase); + if (stripeSparseLen > 0) { + dst.int64().writeArrayBlockWithOffset(blockRowsPath, stripeRows, (int) stripeSparseLen, stripeSparseOffset); + dst.int64().writeArrayBlockWithOffset(blockColsPath, stripeCols, (int) stripeSparseLen, stripeSparseOffset); + dst.float64().writeArrayBlockWithOffset(blockValsPath, stripeVals, (int) stripeSparseLen, stripeSparseOffset); + } + } else { + final var typedBlocks = (LongStripeBlocks) blocks; + final var sparseVals = typedBlocks.sparseVals(); + final var denseFlats = typedBlocks.denseFlats(); + final long[] stripeVals = stripeSparseLen > 0 ? new long[(int) stripeSparseLen] : new long[0]; + for (int i = 0; i < colStripes.length; i++) { + final int blockLen = lengths[i]; + if (blockLen <= 0) { + continue; + } + final int colStripe = colStripes[i]; + + if (denseFlags[i]) { + offsetRow[colStripe] = -denseOffset - 1L; + lengthRow[colStripe] = blockLen; + dst.int64().writeMDArrayBlockWithOffset( + denseBlocksPath, + new MDLongArray(denseFlats[i], new int[]{1, 1, SUBMATRIX_SIZE, SUBMATRIX_SIZE}), + new long[]{denseOffset, 0L, 0L, 0L} + ); + denseOffset++; + } else { + offsetRow[colStripe] = sparseOffset; + lengthRow[colStripe] = blockLen; + final var blockRows = sparseRows[i]; + final var blockCols = sparseCols[i]; + final var blockVals = sparseVals[i]; + System.arraycopy(blockRows, 0, stripeRows, stripePos, blockLen); + System.arraycopy(blockCols, 0, stripeCols, stripePos, blockLen); + System.arraycopy(blockVals, 0, stripeVals, stripePos, blockLen); + stripePos += blockLen; + sparseOffset += blockLen; + } + } + final long rowBase = (long) rowStripe * stripeCount; + dst.int64().writeArrayBlockWithOffset(blockOffsetPath, offsetRow, stripeCount, rowBase); + dst.int64().writeArrayBlockWithOffset(blockLengthPath, lengthRow, stripeCount, rowBase); + if (stripeSparseLen > 0) { + dst.int64().writeArrayBlockWithOffset(blockRowsPath, stripeRows, (int) stripeSparseLen, stripeSparseOffset); + dst.int64().writeArrayBlockWithOffset(blockColsPath, stripeCols, (int) stripeSparseLen, stripeSparseOffset); + dst.int64().writeArrayBlockWithOffset(blockValsPath, stripeVals, (int) stripeSparseLen, stripeSparseOffset); } } - - final long rowBase = (long) rowStripe * stripeCount; - dst.int64().writeArrayBlockWithOffset(blockOffsetPath, offsetRow, stripeCount, rowBase); - dst.int64().writeArrayBlockWithOffset(blockLengthPath, lengthRow, stripeCount, rowBase); - if (stripeSparseLen > 0) { - dst.int64().writeArrayBlockWithOffset(blockRowsPath, stripeRows, (int) stripeSparseLen, stripeSparseOffset); - dst.int64().writeArrayBlockWithOffset(blockColsPath, stripeCols, (int) stripeSparseLen, stripeSparseOffset); - dst.int64().writeArrayBlockWithOffset(blockValsPath, stripeVals, (int) stripeSparseLen, stripeSparseOffset); - } - return new SaveBlockResult(sparseOffset, denseOffset); } @@ -426,6 +476,7 @@ private static SaveBlockResult saveStripeBlocks( final int stripeCount, final long @NotNull [] allRowsStartIndices, final int stripeWorkers, + final boolean floatingPointSignal, final @NotNull java.util.function.IntConsumer countingProgressReporter ) { if (stripeCount <= 0) { @@ -447,7 +498,7 @@ private static SaveBlockResult saveStripeBlocks( for (int batchStart = 0; batchStart < stripeCount; batchStart += batchSize) { checkInterrupted(); final int end = Math.min(stripeCount, batchStart + batchSize); - final var blocks = readStripeBatch(readerHolder.get(), resolution, batchStart, end, allRowsStartIndices); + final var blocks = readStripeBatch(readerHolder.get(), resolution, batchStart, end, allRowsStartIndices, floatingPointSignal); for (int i = 0; i < blocks.length; i++) { final int stripeIdx = batchStart + i; final var block = blocks[i]; @@ -466,7 +517,7 @@ private static SaveBlockResult saveStripeBlocks( final int end = Math.min(stripeCount, batchStart + batchSize); futures.add(stripeExecutor.submit(() -> { checkInterrupted(); - final var blocks = readStripeBatch(readerHolder.get(), resolution, start, end, allRowsStartIndices); + final var blocks = readStripeBatch(readerHolder.get(), resolution, start, end, allRowsStartIndices, floatingPointSignal); for (int i = 0; i < blocks.length; i++) { final int stripeIdx = start + i; final var block = blocks[i]; @@ -553,7 +604,8 @@ private static SaveBlockResult saveStripeBlocks( final long resolution, final int startStripe, final int endStripe, - final long @NotNull [] allRowsStartIndices + final long @NotNull [] allRowsStartIndices, + final boolean floatingPointSignal ) { final int actualEnd = Math.max(startStripe, endStripe); final int count = actualEnd - startStripe; @@ -568,7 +620,9 @@ private static SaveBlockResult saveStripeBlocks( final int length = (int) (endOffset - startOffset); if (length <= 0) { for (int i = 0; i < count; i++) { - blocks[i] = new PixelBlock(new long[0], new long[0], new long[0]); + blocks[i] = floatingPointSignal + ? new DoublePixelBlock(new long[0], new long[0], new double[0]) + : new LongPixelBlock(new long[0], new long[0], new long[0]); } return blocks; } @@ -576,7 +630,8 @@ private static SaveBlockResult saveStripeBlocks( final var base = "/resolutions/" + resolution + "/pixels/"; final var rows = src.int64().readArrayBlockWithOffset(base + "bin1_id", length, startOffset); final var cols = src.int64().readArrayBlockWithOffset(base + "bin2_id", length, startOffset); - final var vals = src.int64().readArrayBlockWithOffset(base + "count", length, startOffset); + final var doubleVals = floatingPointSignal ? src.float64().readArrayBlockWithOffset(base + "count", length, startOffset) : null; + final var longVals = floatingPointSignal ? null : src.int64().readArrayBlockWithOffset(base + "count", length, startOffset); for (int i = 0; i < count; i++) { final int stripeIdx = startStripe + i; @@ -586,12 +641,20 @@ private static SaveBlockResult saveStripeBlocks( final int stripeEnd = (int) (allRowsStartIndices[stripeEndBin] - startOffset); final int stripeLen = Math.max(0, stripeEnd - stripeStart); if (stripeLen <= 0) { - blocks[i] = new PixelBlock(new long[0], new long[0], new long[0]); + blocks[i] = floatingPointSignal + ? new DoublePixelBlock(new long[0], new long[0], new double[0]) + : new LongPixelBlock(new long[0], new long[0], new long[0]); } else { - blocks[i] = new PixelBlock( + blocks[i] = floatingPointSignal + ? new DoublePixelBlock( Arrays.copyOfRange(rows, stripeStart, stripeEnd), Arrays.copyOfRange(cols, stripeStart, stripeEnd), - Arrays.copyOfRange(vals, stripeStart, stripeEnd) + Arrays.copyOfRange(doubleVals, stripeStart, stripeEnd) + ) + : new LongPixelBlock( + Arrays.copyOfRange(rows, stripeStart, stripeEnd), + Arrays.copyOfRange(cols, stripeStart, stripeEnd), + Arrays.copyOfRange(longVals, stripeStart, stripeEnd) ); } } @@ -608,6 +671,103 @@ private static int resolveBatchSize(final int stripeWorkers) { private static @NotNull StripeBlocks buildStripeBlocks( final @NotNull PixelBlock block, final int stripeCount + ) { + if (block instanceof DoublePixelBlock doubleBlock) { + return buildStripeBlocks(doubleBlock, stripeCount); + } + if (block instanceof LongPixelBlock longBlock) { + return buildStripeBlocks(longBlock, stripeCount); + } + throw new IllegalStateException("Unsupported pixel block type: " + block.getClass().getName()); + } + + private static @NotNull StripeBlocks buildStripeBlocks( + final @NotNull DoublePixelBlock block, + final int stripeCount + ) { + final var rows = block.rows(); + final var cols = block.cols(); + final var values = block.values(); + final int n = rows.length; + + final int maxTouched = Math.min(stripeCount, n); + final var buffer = acquireCountBuffer(stripeCount, maxTouched); + final int[] counts = buffer.counts(); + final int[] touched = buffer.touched(); + int touchedCount = 0; + + for (int i = 0; i < n; i++) { + final int colStripe = (int) (cols[i] / SUBMATRIX_SIZE); + if (counts[colStripe]++ == 0) { + touched[touchedCount++] = colStripe; + } + } + + Arrays.sort(touched, 0, touchedCount); + final int nonEmptyBlocks = touchedCount; + final int[] blockColStripes = new int[nonEmptyBlocks]; + final boolean[] denseFlags = new boolean[nonEmptyBlocks]; + final int[] blockLengths = new int[nonEmptyBlocks]; + + for (int i = 0; i < touchedCount; i++) { + final int colStripe = touched[i]; + final int count = counts[colStripe]; + counts[colStripe] = 0; + blockColStripes[i] = colStripe; + blockLengths[i] = count; + denseFlags[i] = count >= DENSE_THRESHOLD; + } + + final long[][] sparseRows = new long[nonEmptyBlocks][]; + final long[][] sparseCols = new long[nonEmptyBlocks][]; + final double[][] sparseVals = new double[nonEmptyBlocks][]; + final double[][] denseFlats = new double[nonEmptyBlocks][]; + + for (int i = 0; i < nonEmptyBlocks; i++) { + if (denseFlags[i]) { + denseFlats[i] = new double[SUBMATRIX_SIZE * SUBMATRIX_SIZE]; + } else { + final int len = blockLengths[i]; + sparseRows[i] = new long[len]; + sparseCols[i] = new long[len]; + sparseVals[i] = new double[len]; + } + } + + final int[] colStripeToIndex = new int[stripeCount]; + Arrays.fill(colStripeToIndex, -1); + for (int i = 0; i < nonEmptyBlocks; i++) { + colStripeToIndex[blockColStripes[i]] = i; + } + + final int[] sparsePositions = new int[nonEmptyBlocks]; + for (int i = 0; i < n; i++) { + final int colStripe = (int) (cols[i] / SUBMATRIX_SIZE); + final int idx = colStripeToIndex[colStripe]; + final int intraRow = (int) (rows[i] % SUBMATRIX_SIZE); + final int intraCol = (int) (cols[i] % SUBMATRIX_SIZE); + if (denseFlags[idx]) { + denseFlats[idx][intraRow * SUBMATRIX_SIZE + intraCol] += values[i]; + } else { + final int pos = sparsePositions[idx]++; + sparseRows[idx][pos] = intraRow; + sparseCols[idx][pos] = intraCol; + sparseVals[idx][pos] = values[i]; + } + } + + for (int i = 0; i < nonEmptyBlocks; i++) { + if (!denseFlags[i] && blockLengths[i] > 1) { + sortSparseBlockRowMajor(sparseRows[i], sparseCols[i], sparseVals[i]); + } + } + + return new DoubleStripeBlocks(blockColStripes, blockLengths, denseFlags, sparseRows, sparseCols, sparseVals, denseFlats); + } + + private static @NotNull StripeBlocks buildStripeBlocks( + final @NotNull LongPixelBlock block, + final int stripeCount ) { final var rows = block.rows(); final var cols = block.cols(); @@ -686,7 +846,45 @@ private static int resolveBatchSize(final int stripeWorkers) { } } - return new StripeBlocks(blockColStripes, blockLengths, denseFlags, sparseRows, sparseCols, sparseVals, denseFlats); + return new LongStripeBlocks(blockColStripes, blockLengths, denseFlags, sparseRows, sparseCols, sparseVals, denseFlats); + } + + private static void sortSparseBlockRowMajor( + final long @NotNull [] rows, + final long @NotNull [] cols, + final double @NotNull [] vals + ) { + final int n = rows.length; + if (n <= 1) { + return; + } + final int bucketSize = SUBMATRIX_SIZE * SUBMATRIX_SIZE; + final int[] counts = new int[bucketSize]; + final int[] keys = new int[n]; + for (int i = 0; i < n; i++) { + final int key = (int) (rows[i] * SUBMATRIX_SIZE + cols[i]); + keys[i] = key; + counts[key]++; + } + int sum = 0; + for (int i = 0; i < bucketSize; i++) { + final int c = counts[i]; + counts[i] = sum; + sum += c; + } + final long[] sortedRows = new long[n]; + final long[] sortedCols = new long[n]; + final double[] sortedVals = new double[n]; + for (int i = 0; i < n; i++) { + final int key = keys[i]; + final int pos = counts[key]++; + sortedRows[pos] = rows[i]; + sortedCols[pos] = cols[i]; + sortedVals[pos] = vals[i]; + } + System.arraycopy(sortedRows, 0, rows, 0, n); + System.arraycopy(sortedCols, 0, cols, 0, n); + System.arraycopy(sortedVals, 0, vals, 0, n); } private static void sortSparseBlockRowMajor( @@ -751,7 +949,9 @@ private static void sortSparseBlockRowMajor( checkInterrupted(); final var block = blocks[i]; if (block.length() > 0) { - sortedBatch[i] = sortStripePixels(block.rows(), block.cols(), block.values()); + sortedBatch[i] = block instanceof DoublePixelBlock doubleBlock + ? sortStripePixels(doubleBlock.rows(), doubleBlock.cols(), doubleBlock.values()) + : sortStripePixels(((LongPixelBlock) block).rows(), ((LongPixelBlock) block).cols(), toDoubleArray(((LongPixelBlock) block).values())); } } return sortedBatch; @@ -764,7 +964,9 @@ private static void sortSparseBlockRowMajor( checkInterrupted(); final var block = blocks[idx]; if (block.length() > 0) { - sortedBatch[idx] = sortStripePixels(block.rows(), block.cols(), block.values()); + sortedBatch[idx] = block instanceof DoublePixelBlock doubleBlock + ? sortStripePixels(doubleBlock.rows(), doubleBlock.cols(), doubleBlock.values()) + : sortStripePixels(((LongPixelBlock) block).rows(), ((LongPixelBlock) block).cols(), toDoubleArray(((LongPixelBlock) block).values())); } })); } @@ -827,22 +1029,26 @@ private static void checkInterrupted() { final int nextBin = Math.min((rowStripe + 1) * SUBMATRIX_SIZE, allRowsStartIndices.length - 1); final long endOffset = allRowsStartIndices[nextBin]; final int length = (int) (endOffset - startOffset); + final var base = "/resolutions/" + resolution + "/pixels/"; + final var floatingPointSignal = isFloatingPointDataset(src, base + "count"); if (length <= 0) { - return new PixelBlock(new long[0], new long[0], new long[0]); + return floatingPointSignal + ? new DoublePixelBlock(new long[0], new long[0], new double[0]) + : new LongPixelBlock(new long[0], new long[0], new long[0]); } - final var base = "/resolutions/" + resolution + "/pixels/"; final var rows = src.int64().readArrayBlockWithOffset(base + "bin1_id", length, startOffset); final var cols = src.int64().readArrayBlockWithOffset(base + "bin2_id", length, startOffset); - final var vals = src.int64().readArrayBlockWithOffset(base + "count", length, startOffset); - return new PixelBlock(rows, cols, vals); + return floatingPointSignal + ? new DoublePixelBlock(rows, cols, src.float64().readArrayBlockWithOffset(base + "count", length, startOffset)) + : new LongPixelBlock(rows, cols, src.int64().readArrayBlockWithOffset(base + "count", length, startOffset)); } private static @NotNull SortedStripePixels sortStripePixels( final long @NotNull [] rows, final long @NotNull [] cols, - final long @NotNull [] values + final double @NotNull [] values ) { final int n = rows.length; final var order = new Integer[n]; @@ -860,7 +1066,7 @@ private static void checkInterrupted() { final var sortedColStripes = new long[n]; final var sortedIntraRows = new int[n]; final var sortedIntraCols = new int[n]; - final var sortedValues = new long[n]; + final var sortedValues = new double[n]; for (int i = 0; i < n; i++) { final int srcIdx = order[i]; @@ -1163,6 +1369,74 @@ private static void runParallelFor(final int parallelism, final int itemCount, f } } + private static boolean isFloatingPointDataset(final @NotNull IHDF5Reader reader, final @NotNull String path) { + return reader.object().getDataSetInformation(path).getTypeInformation().getDataClass() == HDF5DataClass.FLOAT; + } + + private static double @NotNull [] readNumericArrayBlockWithOffset( + final @NotNull IHDF5Reader reader, + final @NotNull String path, + final int length, + final long offset + ) { + if (isFloatingPointDataset(reader, path)) { + return reader.float64().readArrayBlockWithOffset(path, length, offset); + } + return toDoubleArray(reader.int64().readArrayBlockWithOffset(path, length, offset)); + } + + private static void createNumericArray( + final @NotNull IHDF5Writer writer, + final @NotNull String path, + final long length, + final int chunkSize, + final boolean floatingPoint, + final @NotNull HDF5IntStorageFeatures intStorageFeatures, + final @NotNull HDF5FloatStorageFeatures floatStorageFeatures + ) { + if (floatingPoint) { + writer.float64().createArray(path, length, safeChunkLen(length, chunkSize), floatStorageFeatures); + } else { + writer.int64().createArray(path, length, safeChunkLen(length, chunkSize), intStorageFeatures); + } + } + + private static void createNumericMDArray( + final @NotNull IHDF5Writer writer, + final @NotNull String path, + final long @NotNull [] dimensions, + final int @NotNull [] blockDimensions, + final boolean floatingPoint, + final @NotNull HDF5IntStorageFeatures intStorageFeatures, + final @NotNull HDF5FloatStorageFeatures floatStorageFeatures + ) { + if (floatingPoint) { + writer.float64().createMDArray(path, dimensions, blockDimensions, floatStorageFeatures); + } else { + writer.int64().createMDArray(path, dimensions, blockDimensions, intStorageFeatures); + } + } + + private static double @NotNull [] toDoubleArray(final long @NotNull [] input) { + final var result = new double[input.length]; + for (int i = 0; i < input.length; i++) { + result[i] = input[i]; + } + return result; + } + + private static long @NotNull [] toLongArray(final double @NotNull [] input) { + final var result = new long[input.length]; + for (int i = 0; i < input.length; i++) { + result[i] = Math.round(input[i]); + } + return result; + } + + private static @NotNull MDLongArray toLongMdArray(final double @NotNull [] input) { + return new MDLongArray(toLongArray(input), new int[]{1, 1, SUBMATRIX_SIZE, SUBMATRIX_SIZE}); + } + private static @NotNull String resolveNameLengthPath(final @NotNull IHDF5Reader src, final long resolution) { final var perResolution = "/resolutions/" + resolution + "/chroms"; if (src.object().isGroup(perResolution) && src.object().isDataSet(perResolution + "/name") && src.object().isDataSet(perResolution + "/length")) { @@ -1290,10 +1564,26 @@ private void finish() { } private static final StripeBlocks EMPTY_BLOCKS = - new StripeBlocks(new int[0], new int[0], new boolean[0], new long[0][], new long[0][], new long[0][], new long[0][]); + new LongStripeBlocks(new int[0], new int[0], new boolean[0], new long[0][], new long[0][], new long[0][], new long[0][]); + + private sealed interface PixelBlock permits LongPixelBlock, DoublePixelBlock { + long @NotNull [] rows(); + + long @NotNull [] cols(); + + int length(); + } - private record PixelBlock(long @NotNull [] rows, long @NotNull [] cols, long @NotNull [] values) { - int length() { + private record LongPixelBlock(long @NotNull [] rows, long @NotNull [] cols, long @NotNull [] values) implements PixelBlock { + @Override + public int length() { + return rows.length; + } + } + + private record DoublePixelBlock(long @NotNull [] rows, long @NotNull [] cols, double @NotNull [] values) implements PixelBlock { + @Override + public int length() { return rows.length; } } @@ -1302,11 +1592,23 @@ private record SortedStripePixels( long @NotNull [] colStripes, int @NotNull [] intraRows, int @NotNull [] intraCols, - long @NotNull [] values + double @NotNull [] values ) { } - private record StripeBlocks( + private sealed interface StripeBlocks permits LongStripeBlocks, DoubleStripeBlocks { + int @NotNull [] colStripes(); + + int @NotNull [] blockLengths(); + + boolean @NotNull [] denseFlags(); + + long @NotNull [] @NotNull [] sparseRows(); + + long @NotNull [] @NotNull [] sparseCols(); + } + + private record LongStripeBlocks( int @NotNull [] colStripes, int @NotNull [] blockLengths, boolean @NotNull [] denseFlags, @@ -1314,7 +1616,18 @@ private record StripeBlocks( long @NotNull [] @NotNull [] sparseCols, long @NotNull [] @NotNull [] sparseVals, long @NotNull [] @NotNull [] denseFlats - ) { + ) implements StripeBlocks { + } + + private record DoubleStripeBlocks( + int @NotNull [] colStripes, + int @NotNull [] blockLengths, + boolean @NotNull [] denseFlags, + long @NotNull [] @NotNull [] sparseRows, + long @NotNull [] @NotNull [] sparseCols, + double @NotNull [] @NotNull [] sparseVals, + double @NotNull [] @NotNull [] denseFlats + ) implements StripeBlocks { } private record StripeCounts(long sparseElementCount, long denseBlockCount) { diff --git a/src/main/java/ru/itmo/ctlab/hict/hict_library/visualization/RawTileWithWeights.java b/src/main/java/ru/itmo/ctlab/hict/hict_library/visualization/RawTileWithWeights.java index 6e026c6..8cf5069 100644 --- a/src/main/java/ru/itmo/ctlab/hict/hict_library/visualization/RawTileWithWeights.java +++ b/src/main/java/ru/itmo/ctlab/hict/hict_library/visualization/RawTileWithWeights.java @@ -27,6 +27,6 @@ import org.jetbrains.annotations.NotNull; import org.jetbrains.annotations.Nullable; -public record RawTileWithWeights(long @NotNull [][] values, double @Nullable [] rowWeights, +public record RawTileWithWeights(double @NotNull [][] values, double @Nullable [] rowWeights, double @Nullable [] columnWeights) { } diff --git a/src/main/java/ru/itmo/ctlab/hict/hict_library/visualization/TileVisualizationProcessor.java b/src/main/java/ru/itmo/ctlab/hict/hict_library/visualization/TileVisualizationProcessor.java index 2f3f81d..a329532 100644 --- a/src/main/java/ru/itmo/ctlab/hict/hict_library/visualization/TileVisualizationProcessor.java +++ b/src/main/java/ru/itmo/ctlab/hict/hict_library/visualization/TileVisualizationProcessor.java @@ -37,7 +37,6 @@ import java.util.concurrent.Executors; import java.util.concurrent.atomic.AtomicInteger; import java.util.stream.DoubleStream; -import java.util.stream.LongStream; import java.util.stream.Stream; @RequiredArgsConstructor @@ -58,21 +57,6 @@ protected DoubleStream applyCoolerWeightsToRow(final DoubleStream rowStream, fin return doubleStream; } - protected DoubleStream applyCoolerWeightsToRow(final LongStream rowStream, final double rowWeight, final double @Nullable [] columnWeights) { - DoubleStream doubleStream; - if (rowWeight != 1.0) { - doubleStream = rowStream.mapToDouble(signal -> signal * rowWeight); - } else { - doubleStream = rowStream.mapToDouble(signal -> (double) signal); - } - - if (columnWeights != null) { - final var atomicColumnIndex = new AtomicInteger(); - doubleStream = doubleStream.sequential().map(signal -> signal * columnWeights[atomicColumnIndex.getAndIncrement()]); - } - return doubleStream; - } - public @NotNull TileWithWeights applyWeightsToTile(final @NotNull RawTileWithWeights rawTile) { final var rowCount = rawTile.values().length; final var columnCount = (rawTile.values().length > 0) ? (rawTile.values()[0].length) : 0; @@ -99,7 +83,7 @@ protected DoubleStream applyCoolerWeightsToRow(final LongStream rowStream, final final var rawRow = rawValues[rowIndex]; final var wRow = new double[columnCount]; for (int j = 0; j < columnCount; j++) { - wRow[j] = ((double) rawRow[j]) * rowWeights[rowIndex] * colWeights[j]; + wRow[j] = rawRow[j] * rowWeights[rowIndex] * colWeights[j]; } weightedValues[rowIndex] = wRow; }); @@ -116,22 +100,29 @@ public TileWithWeights processTile(final @NotNull MatrixQueries.MatrixWithWeight final var input = rawTile.matrix(); final var rowWeights = rawTile.rowWeights(); final var columnWeights = rawTile.colWeights(); - final var rowCount = input.length; - final var columnCount = (rowCount > 0) ? input[0].length : 0; + final var rowCount = input.rows(); + final var columnCount = input.cols(); final var result = new double[rowCount][columnCount]; final var resolutionScalingCoeffs = this.chunkedFile.getResolutionScalingCoefficient(); final var resolutionLinearScalingCoeffs = this.chunkedFile.getResolutionLinearScalingCoefficient(); final var resolutionScalingCoeff = resolutionScalingCoeffs[rawTile.resolutionDescriptor().getResolutionOrderInArray()]; final var resolutionLinearScalingCoeff = resolutionLinearScalingCoeffs[rawTile.resolutionDescriptor().getResolutionOrderInArray()]; - for (int rowIndex = 0; rowIndex < input.length; ++rowIndex) { - final var startStream = Arrays.stream(input[rowIndex]).parallel(); + for (int rowIndex = 0; rowIndex < rowCount; ++rowIndex) { + final DoubleStream startStream; + if (input instanceof MatrixQueries.DoubleMatrix doubleMatrix) { + startStream = Arrays.stream(doubleMatrix.values()[rowIndex]).parallel(); + } else if (input instanceof MatrixQueries.LongMatrix longMatrix) { + startStream = Arrays.stream(longMatrix.values()[rowIndex]).parallel().mapToDouble(value -> value); + } else { + throw new IllegalStateException("Unsupported raw matrix type: " + input.getClass().getName()); + } final var rowWeight = (rowWeights != null) ? rowWeights[rowIndex] : 1.0; final var pre = visualizationOptions.getLnPreLogBase(); final var post = visualizationOptions.getLnPostLogBase(); - DoubleStream doubleStream = startStream.mapToDouble(signal -> (double) signal); + DoubleStream doubleStream = startStream; if (pre > 0) { doubleStream = doubleStream.map(signal -> Math.log1p(signal) / pre); @@ -160,8 +151,8 @@ public TileWithWeights processTile(final @NotNull MatrixQueries.MatrixWithWeight public @NotNull BufferedImage visualizeTile(final @NotNull MatrixQueries.MatrixWithWeights rawTile, final @NotNull SimpleVisualizationOptions options) { final var input = rawTile.matrix(); - final var rowCount = input.length; - final var columnCount = (rowCount > 0) ? input[0].length : 0; + final var rowCount = input.rows(); + final var columnCount = input.cols(); final var normalized = processTile(rawTile, options); final var colormap = options.getColormap(); final var boxedARGBValues = Arrays.stream(normalized.values()) diff --git a/src/main/java/ru/itmo/ctlab/hict/hict_server/handlers/tiles/TileHandlersHolder.java b/src/main/java/ru/itmo/ctlab/hict/hict_server/handlers/tiles/TileHandlersHolder.java index e9f4846..cfa68d8 100644 --- a/src/main/java/ru/itmo/ctlab/hict/hict_server/handlers/tiles/TileHandlersHolder.java +++ b/src/main/java/ru/itmo/ctlab/hict/hict_server/handlers/tiles/TileHandlersHolder.java @@ -574,8 +574,9 @@ private MatrixResponsePayload computeMatrixQueryResponse(final @NotNull JsonObje endColPx, true ); - final var rowCount = matrixWithWeights.matrix().length; - final var columnCount = rowCount > 0 ? matrixWithWeights.matrix()[0].length : 0; + final var rawMatrix = matrixWithWeights.matrix(); + final var rowCount = rawMatrix.rows(); + final var columnCount = rawMatrix.cols(); final long elementCount = (long) rowCount * (long) columnCount; if (elementCount > MAX_MATRIX_QUERY_ELEMENTS) { throw new IllegalArgumentException( @@ -588,7 +589,6 @@ private MatrixResponsePayload computeMatrixQueryResponse(final @NotNull JsonObje final var includeWeights = request.getBoolean("includeWeights", false); final double[][] signalMatrix; - final long[][] rawMatrix = matrixWithWeights.matrix(); switch (signalMode) { case RAW_COUNTS -> signalMatrix = null; case COOLER_WEIGHTED -> signalMatrix = computeCoolerWeightedSignal(rawMatrix, matrixWithWeights.rowWeights(), matrixWithWeights.colWeights()); @@ -636,8 +636,15 @@ private MatrixResponsePayload computeMatrixQueryResponse(final @NotNull JsonObje .put("startColPx", matrixWithWeights.startColIncl()) .put("endColPx", matrixWithWeights.startColIncl() + columnCount); if (signalMode == MatrixSignalMode.RAW_COUNTS) { - payload.put("dtype", "int64"); - payload.put("values", flattenLongMatrix(rawMatrix, rowCount, columnCount)); + if (rawMatrix instanceof ru.itmo.ctlab.hict.hict_library.chunkedfile.MatrixQueries.DoubleMatrix doubleMatrix) { + payload.put("dtype", "float64"); + payload.put("values", flattenDoubleMatrix(doubleMatrix.values(), rowCount, columnCount)); + } else if (rawMatrix instanceof ru.itmo.ctlab.hict.hict_library.chunkedfile.MatrixQueries.LongMatrix longMatrix) { + payload.put("dtype", "int64"); + payload.put("values", flattenLongMatrix(longMatrix.values(), rowCount, columnCount)); + } else { + throw new IllegalStateException("Unsupported raw matrix type: " + rawMatrix.getClass().getName()); + } } else { payload.put("dtype", "float64"); payload.put("values", flattenDoubleMatrix(signalMatrix, rowCount, columnCount)); @@ -652,23 +659,29 @@ private MatrixResponsePayload computeMatrixQueryResponse(final @NotNull JsonObje switch (format) { case BINARY_INT64 -> { - headers.put("x-hict-dtype", "int64"); - final var binary = encodeLongMatrixLittleEndian(rawMatrix, rowCount, columnCount); + final byte[] binary; + if (signalMode == MatrixSignalMode.RAW_COUNTS) { + if (rawMatrix instanceof ru.itmo.ctlab.hict.hict_library.chunkedfile.MatrixQueries.LongMatrix longMatrix) { + headers.put("x-hict-dtype", "int64"); + binary = encodeLongMatrixLittleEndian(longMatrix.values(), rowCount, columnCount); + } else { + throw new IllegalArgumentException("BINARY_INT64 is not supported for float-backed RAW_COUNTS matrices"); + } + } else { + headers.put("x-hict-dtype", "int64"); + binary = encodeLongMatrixLittleEndian(signalMatrix, rowCount, columnCount); + } return new MatrixResponsePayload("application/octet-stream", null, Buffer.buffer(binary), headers); } case BINARY_FLOAT64 -> { headers.put("x-hict-dtype", "float64"); - final var source = signalMode == MatrixSignalMode.RAW_COUNTS - ? asDoubleMatrix(rawMatrix, rowCount, columnCount) - : signalMatrix; + final var source = signalMode == MatrixSignalMode.RAW_COUNTS ? toDoubleMatrix(rawMatrix) : signalMatrix; final var binary = encodeDoubleMatrixLittleEndian(source, rowCount, columnCount); return new MatrixResponsePayload("application/octet-stream", null, Buffer.buffer(binary), headers); } case BINARY_FLOAT32 -> { headers.put("x-hict-dtype", "float32"); - final var source = signalMode == MatrixSignalMode.RAW_COUNTS - ? asDoubleMatrix(rawMatrix, rowCount, columnCount) - : signalMatrix; + final var source = signalMode == MatrixSignalMode.RAW_COUNTS ? toDoubleMatrix(rawMatrix) : signalMatrix; final var binary = encodeFloatMatrixLittleEndian(source, rowCount, columnCount); return new MatrixResponsePayload("application/octet-stream", null, Buffer.buffer(binary), headers); } @@ -676,17 +689,17 @@ private MatrixResponsePayload computeMatrixQueryResponse(final @NotNull JsonObje } } - private static double[][] computeCoolerWeightedSignal(final long[][] rawMatrix, + private static double[][] computeCoolerWeightedSignal(final ru.itmo.ctlab.hict.hict_library.chunkedfile.MatrixQueries.RawMatrix rawMatrix, final double[] rowWeights, final double[] colWeights) { - final var rowCount = rawMatrix.length; - final var columnCount = rowCount > 0 ? rawMatrix[0].length : 0; + final var rowCount = rawMatrix.rows(); + final var columnCount = rawMatrix.cols(); final var result = new double[rowCount][columnCount]; for (int row = 0; row < rowCount; row++) { final var rowWeight = rowWeights != null && row < rowWeights.length ? rowWeights[row] : 1.0d; for (int col = 0; col < columnCount; col++) { final var colWeight = colWeights != null && col < colWeights.length ? colWeights[col] : 1.0d; - result[row][col] = rawMatrix[row][col] * rowWeight * colWeight; + result[row][col] = rawMatrix.getAsDouble(row, col) * rowWeight * colWeight; } } return result; @@ -700,8 +713,8 @@ private double[][] computePipelineSignalMatrix(final @NotNull ChunkedFile primar final @NotNull RenderPipelineConfig pipelineConfig, final Track1DManager track1DManager) { final var primaryValues = primaryMatrixWithWeights.matrix(); - final var rowCount = primaryValues.length; - final var columnCount = rowCount > 0 ? primaryValues[0].length : 0; + final var rowCount = primaryValues.rows(); + final var columnCount = primaryValues.cols(); final var result = new double[rowCount][columnCount]; if (rowCount == 0 || columnCount == 0) { return result; @@ -710,9 +723,8 @@ private double[][] computePipelineSignalMatrix(final @NotNull ChunkedFile primar final var secondaryValues = new double[rowCount][columnCount]; if (secondaryMatrixWithWeights != null && secondaryChunkedFile != null) { final var candidate = secondaryMatrixWithWeights.matrix(); - final var candidateRowCount = candidate.length; - final var candidateColCount = - candidateRowCount > 0 && candidate[0] != null ? candidate[0].length : 0; + final var candidateRowCount = candidate.rows(); + final var candidateColCount = candidate.cols(); final var rowOffset = (int) (secondaryMatrixWithWeights.startRowIncl() - primaryMatrixWithWeights.startRowIncl()); final var colOffset = @@ -722,17 +734,12 @@ private double[][] computePipelineSignalMatrix(final @NotNull ChunkedFile primar if (dstRow < 0 || dstRow >= rowCount) { continue; } - final var sourceRow = candidate[row]; - if (sourceRow == null) { - continue; - } - final var sourceColCount = Math.min(sourceRow.length, candidateColCount); - for (int col = 0; col < sourceColCount; col++) { + for (int col = 0; col < candidateColCount; col++) { final var dstCol = col + colOffset; if (dstCol < 0 || dstCol >= columnCount) { continue; } - secondaryValues[dstRow][dstCol] = sourceRow[col]; + secondaryValues[dstRow][dstCol] = candidate.getAsDouble(row, col); } } } @@ -849,7 +856,7 @@ private double[][] computePipelineSignalMatrix(final @NotNull ChunkedFile primar final var rowBin = rowBinValues[row]; final var rowBp = rowBpValues[row]; for (int col = 0; col < columnCount; ++col) { - final var primaryValue = (double) primaryValues[row][col]; + final var primaryValue = primaryValues.getAsDouble(row, col); final var secondaryValue = secondaryValues[row][col]; final var colPx = colPxValues[col]; final var colBin = colBinValues[col]; @@ -916,16 +923,15 @@ private static byte[] encodeLongMatrixLittleEndian(final long[][] matrix, final return bb.array(); } - private static double[][] asDoubleMatrix(final long[][] source, final int rows, final int cols) { - final var result = new double[rows][cols]; + private static byte[] encodeLongMatrixLittleEndian(final double[][] matrix, final int rows, final int cols) { + final var bb = ByteBuffer.allocate(rows * cols * Long.BYTES).order(ByteOrder.LITTLE_ENDIAN); for (int row = 0; row < rows; row++) { - final var src = source[row]; - final var dst = result[row]; + final var sourceRow = matrix[row]; for (int col = 0; col < cols; col++) { - dst[col] = src[col]; + bb.putLong(Math.round(sourceRow[col])); } } - return result; + return bb.array(); } private static ArrayList toJsonArray(final double[] values, final int expectedLength) { @@ -937,8 +943,8 @@ private static ArrayList toJsonArray(final double[] values, final int ex return result; } - private static ArrayList flattenLongMatrix(final long[][] matrix, final int rows, final int cols) { - final var result = new ArrayList(rows * cols); + private static ArrayList flattenDoubleMatrix(final double[][] matrix, final int rows, final int cols) { + final var result = new ArrayList(rows * cols); for (int row = 0; row < rows; row++) { final var sourceRow = matrix[row]; for (int col = 0; col < cols; col++) { @@ -948,8 +954,8 @@ private static ArrayList flattenLongMatrix(final long[][] matrix, final in return result; } - private static ArrayList flattenDoubleMatrix(final double[][] matrix, final int rows, final int cols) { - final var result = new ArrayList(rows * cols); + private static ArrayList flattenLongMatrix(final long[][] matrix, final int rows, final int cols) { + final var result = new ArrayList(rows * cols); for (int row = 0; row < rows; row++) { final var sourceRow = matrix[row]; for (int col = 0; col < cols; col++) { @@ -959,6 +965,27 @@ private static ArrayList flattenDoubleMatrix(final double[][] matrix, fi return result; } + private static double[][] toDoubleMatrix(final ru.itmo.ctlab.hict.hict_library.chunkedfile.MatrixQueries.RawMatrix matrix) { + if (matrix instanceof ru.itmo.ctlab.hict.hict_library.chunkedfile.MatrixQueries.DoubleMatrix doubleMatrix) { + return doubleMatrix.values(); + } + if (matrix instanceof ru.itmo.ctlab.hict.hict_library.chunkedfile.MatrixQueries.LongMatrix longMatrix) { + final var rows = longMatrix.rows(); + final var cols = longMatrix.cols(); + final var result = new double[rows][cols]; + final var source = longMatrix.values(); + for (int row = 0; row < rows; row++) { + final var sourceRow = source[row]; + final var targetRow = result[row]; + for (int col = 0; col < cols; col++) { + targetRow[col] = sourceRow[col]; + } + } + return result; + } + throw new IllegalStateException("Unsupported raw matrix type: " + matrix.getClass().getName()); + } + private long resolveRangeStart(final @NotNull JsonObject request, final @NotNull Axis axis, final @NotNull QueryLengthUnit units) { @@ -1210,8 +1237,8 @@ public enum MatrixResponseFormat { final @NotNull RenderPipelineConfig pipelineConfig, final Track1DManager track1DManager) { final var primaryValues = primaryMatrixWithWeights.matrix(); - final var rowCount = primaryValues.length; - final var columnCount = rowCount > 0 ? primaryValues[0].length : 0; + final var rowCount = primaryValues.rows(); + final var columnCount = primaryValues.cols(); final var image = new BufferedImage(columnCount, rowCount, BufferedImage.TYPE_INT_ARGB); final var rgba = new int[Math.max(0, rowCount * columnCount)]; if (rowCount == 0 || columnCount == 0) { @@ -1221,9 +1248,8 @@ public enum MatrixResponseFormat { final var secondaryValues = new double[rowCount][columnCount]; if (secondaryMatrixWithWeights != null && secondaryChunkedFile != null) { final var candidate = secondaryMatrixWithWeights.matrix(); - final var candidateRowCount = candidate.length; - final var candidateColCount = - candidateRowCount > 0 && candidate[0] != null ? candidate[0].length : 0; + final var candidateRowCount = candidate.rows(); + final var candidateColCount = candidate.cols(); final var rowOffset = (int) (secondaryMatrixWithWeights.startRowIncl() - primaryMatrixWithWeights.startRowIncl()); final var colOffset = @@ -1233,17 +1259,12 @@ public enum MatrixResponseFormat { if (dstRow < 0 || dstRow >= rowCount) { continue; } - final var sourceRow = candidate[row]; - if (sourceRow == null) { - continue; - } - final var sourceColCount = Math.min(sourceRow.length, candidateColCount); - for (int col = 0; col < sourceColCount; col++) { + for (int col = 0; col < candidateColCount; col++) { final var dstCol = col + colOffset; if (dstCol < 0 || dstCol >= columnCount) { continue; } - secondaryValues[dstRow][dstCol] = sourceRow[col]; + secondaryValues[dstRow][dstCol] = candidate.getAsDouble(row, col); } } } @@ -1361,7 +1382,7 @@ public enum MatrixResponseFormat { final var rowBin = rowBinValues[row]; final var rowBp = rowBpValues[row]; for (int col = 0; col < columnCount; ++col) { - final var primaryValue = (double) primaryValues[row][col]; + final var primaryValue = primaryValues.getAsDouble(row, col); final var secondaryValue = secondaryValues[row][col]; final var colPx = colPxValues[col]; final var colBin = colBinValues[col]; diff --git a/version.txt b/version.txt index ce4bd36..4fececb 100644 --- a/version.txt +++ b/version.txt @@ -1 +1 @@ -1.0.122-11a87a5-webui_936867c \ No newline at end of file +1.0.124-bbee3e5-webui_92b388b \ No newline at end of file