us1zcljz
/
mathnet-numerics
forked from tsai/mathnet-numerics
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity
Browse Source

Data: Matlab: refactoring, simplification, optimizations

provider
Christoph Ruegg 12 years ago
parent
2b177f935a
commit
db8b64c782
10 changed files with 901 additions and 1168 deletions
Whitespace
Show all changes Ignore whitespace when comparing lines Ignore changes in amount of whitespace Ignore changes in whitespace at EOL
Unified View
Diff Options
Show Stats Download Patch File Download Diff File
  1. 649
      src/Data/Matlab/Formatter.cs
  2. 4
      src/Data/Matlab/Matlab.csproj
  3. 6
      src/Data/Matlab/MatlabMatrix.cs
  4. 2
      src/Data/Matlab/MatlabReader.cs
  5. 2
      src/Data/Matlab/MatlabWriter.cs
  6. 150
      src/Data/Matlab/NumericArrayFormatter.cs
  7. 257
      src/Data/Matlab/NumericArrayReader.cs
  8. 447
      src/Data/Matlab/Parser.cs
  9. 301
      src/Data/Matlab/SparseArrayFormatter.cs
  10. 251
      src/Data/Matlab/SparseArrayReader.cs

649
src/Data/Matlab/Formatter.cs
View File

@ -28,7 +28,6 @@
// OTHER DEALINGS IN THE SOFTWARE. // OTHER DEALINGS IN THE SOFTWARE.
// </copyright> // </copyright>
using System; using System;
using System.Collections.Generic; using System.Collections.Generic;
using System.IO; using System.IO;
@ -51,11 +50,6 @@ namespace MathNet.Numerics.Data.Matlab
/// </summary> /// </summary>
const string HeaderText = "MATLAB 5.0 MAT-file, Platform: .NET 4 - Math.NET Numerics, Created on: "; const string HeaderText = "MATLAB 5.0 MAT-file, Platform: .NET 4 - Math.NET Numerics, Created on: ";
/// <summary>
/// The length of the header text.
/// </summary>
const int HeaderTextLength = 116;
/// <summary> /// <summary>
/// Format a matrix block byte array /// Format a matrix block byte array
/// </summary> /// </summary>
@ -77,490 +71,104 @@ namespace MathNet.Numerics.Data.Matlab
throw new ArgumentException(string.Format(Resources.NameCannotContainASpace, name), "name"); throw new ArgumentException(string.Format(Resources.NameCannotContainASpace, name), "name");
} }
if (typeof(T) == typeof(double)) var typeT = typeof (T);
{ bool sparse = matrix.Storage.GetType().GetGenericTypeDefinition() == typeof (SparseCompressedRowMatrixStorage<>);
var sparse = matrix as LinearAlgebra.Double.SparseMatrix; bool doublePrecision = typeT == typeof (double) || typeT == typeof (Complex);
return sparse != null bool complex = typeT == typeof (Complex) || typeT == typeof (Complex32);
? GetSparseDataArray(sparse, name)
: GetDenseDataArray((LinearAlgebra.Double.Matrix)(object)matrix, name);
}
if (typeof(T) == typeof(float))
{
var sparse = matrix as LinearAlgebra.Single.SparseMatrix;
return sparse != null
? GetSparseDataArray(sparse, name)
: GetDenseDataArray((LinearAlgebra.Single.Matrix)(object)matrix, name);
}
if (typeof(T) == typeof(Complex))
{
var sparse = matrix as LinearAlgebra.Complex.SparseMatrix;
return sparse != null
? GetSparseDataArray(sparse, name)
: GetDenseDataArray((LinearAlgebra.Complex.Matrix)(object)matrix, name);
}
if (typeof(T) == typeof(Complex32))
{
var sparse = matrix as LinearAlgebra.Complex32.SparseMatrix;
return sparse != null
? GetSparseDataArray(sparse, name)
: GetDenseDataArray((LinearAlgebra.Complex32.Matrix)(object)matrix, name);
}
throw new NotSupportedException();
}
/// <summary>
/// Writes all matrix blocks to a stream.
/// </summary>
internal static void FormatFile(Stream stream, IEnumerable<MatlabMatrix> matrices)
{
using (var buffer = new BufferedStream(stream))
using (var writer = new BinaryWriter(buffer))
{
WriteHeader(writer);
foreach (var matrix in matrices)
{
// write data type
writer.Write((int)DataType.Compressed);
WriteCompressedData(writer, matrix.Data);
}
writer.Flush();
writer.Close();
}
}
/// <summary>
/// Writes the matrix tag and name.
/// </summary>
/// <param name="writer">The writer we are using.</param>
/// <param name="arrayClass">The array class we are writing.</param>
/// <param name="isComplex">if set to <c>true</c> if this a complex matrix.</param>
/// <param name="name">The name of the matrix.</param>
/// <param name="rows">The number of rows.</param>
/// <param name="columns">The columns of columns.</param>
/// <param name="nzmax">The maximum number of non-zero elements.</param>
static void WriteMatrixTagAndName(BinaryWriter writer, ArrayClass arrayClass, bool isComplex,
string name, int rows, int columns, int nzmax)
{
writer.Write((int)DataType.Matrix);
// add place holder for data size
writer.Write(0);
// write flag, data type and size
writer.Write((int)DataType.UInt32);
writer.Write(8);
// write array class and flags
writer.Write((byte)arrayClass);
if (isComplex)
{
writer.Write((byte)ArrayFlags.Complex);
}
else
{
writer.Write((byte)0);
}
writer.Write((short)0);
writer.Write(nzmax);
// write dimensions
writer.Write((int)DataType.Int32);
writer.Write(8);
writer.Write(rows);
writer.Write(columns);
var nameBytes = Encoding.ASCII.GetBytes(name);
// write name
if (nameBytes.Length > 4)
{
writer.Write((int)DataType.Int8);
writer.Write(nameBytes.Length);
writer.Write(nameBytes);
var pad = 8 - (nameBytes.Length%8);
PadData(writer, pad);
}
else
{
writer.Write((short)DataType.Int8);
writer.Write((short)nameBytes.Length);
writer.Write(nameBytes);
PadData(writer, 4 - nameBytes.Length);
}
}
/// <summary>
/// Gets the dense data array.
/// </summary>
/// <param name="matrix">The matrix to get the data from.</param>
/// <param name="name">The name of the matrix.</param>
/// <returns>The matrix data as an array.</returns>
static MatlabMatrix GetDenseDataArray(Matrix<double> matrix, string name)
{
using (var stream = new MemoryStream())
using (var writer = new BinaryWriter(stream))
{
WriteMatrixTagAndName(writer, ArrayClass.Double, false, name, matrix.RowCount, matrix.ColumnCount, 0);
// write data
writer.Write((int)DataType.Double);
writer.Write(matrix.RowCount*matrix.ColumnCount*8);
for (var j = 0; j < matrix.ColumnCount; j++)
{
var column = matrix.Column(j);
foreach (var value in column)
{
writer.Write(value);
}
}
writer.Flush();
return new MatlabMatrix(name, stream.ToArray());
}
}
/// <summary> int sparseNonZeroValues = 0;
/// Gets the dense data array. if (sparse)
/// </summary>
/// <param name="matrix">The matrix to get the data from.</param>
/// <param name="name">The name of the matrix.</param>
/// <returns>The matrix data as an array.</returns>
static MatlabMatrix GetDenseDataArray(Matrix<float> matrix, string name)
{
using (var stream = new MemoryStream())
using (var writer = new BinaryWriter(stream))
{ {
WriteMatrixTagAndName(writer, ArrayClass.Single, false, name, matrix.RowCount, matrix.ColumnCount, 0); var sparseStorage = (SparseCompressedRowMatrixStorage<T>)matrix.Storage;
sparseNonZeroValues = sparseStorage.ValueCount;
// write data
int size = matrix.RowCount*matrix.ColumnCount*4;
writer.Write((int)DataType.Single);
writer.Write(size);
for (var j = 0; j < matrix.ColumnCount; j++)
{
var column = matrix.Column(j);
foreach (var value in column)
{
writer.Write(value);
}
}
PadData(writer, size%8);
writer.Flush();
return new MatlabMatrix(name, stream.ToArray());
} }
}
/// <summary>
/// Gets the dense data array.
/// </summary>
/// <param name="matrix">The matrix to get the data from.</param>
/// <param name="name">The name of the matrix.</param>
/// <returns>The matrix data as an array.</returns>
static MatlabMatrix GetDenseDataArray(Matrix<Complex> matrix, string name)
{
using (var stream = new MemoryStream()) using (var stream = new MemoryStream())
using (var writer = new BinaryWriter(stream)) using (var writer = new BinaryWriter(stream))
{ {
WriteMatrixTagAndName(writer, ArrayClass.Double, true, name, matrix.RowCount, matrix.ColumnCount, 0); // Array Flags tag: data type + size (8 bytes)
writer.Write((int)DataType.UInt32);
writer.Write(8);
// write data // Array Flags data: flags (byte 3), class (byte 4) (8 bytes)
int size = matrix.RowCount*matrix.ColumnCount*8; writer.Write((byte)(sparse ? ArrayClass.Sparse : doublePrecision ? ArrayClass.Double : ArrayClass.Single));
writer.Write((int)DataType.Double); writer.Write((byte)(complex ? ArrayFlags.Complex : 0));
writer.Write(size); writer.Write((short)0);
writer.Write((int)sparseNonZeroValues);
for (var j = 0; j < matrix.ColumnCount; j++) // Dimensions Array tag: data type + size (8 bytes)
{ writer.Write((int)DataType.Int32);
var column = matrix.Column(j); writer.Write(8);
foreach (var value in column)
{
writer.Write(value.Real);
}
}
writer.Write((int)DataType.Double); // Dimensions Array data: row and column count (8 bytes)
writer.Write(size); writer.Write(matrix.RowCount);
writer.Write(matrix.ColumnCount);
for (var j = 0; j < matrix.ColumnCount; j++) // Array Name:
var nameBytes = Encoding.ASCII.GetBytes(name);
if (nameBytes.Length > 4)
{ {
var column = matrix.Column(j); // long format
foreach (var value in column) writer.Write((int)DataType.Int8);
{ writer.Write(nameBytes.Length);
writer.Write(value.Imaginary); writer.Write(nameBytes);
} PadData(writer, 8 - (nameBytes.Length%8));
} }
else
writer.Flush();
return new MatlabMatrix(name, stream.ToArray());
}
}
/// <summary>
/// Gets the dense data array.
/// </summary>
/// <param name="matrix">The matrix to get the data from.</param>
/// <param name="name">The name of the matrix.</param>
/// <returns>The matrix data as an array.</returns>
static MatlabMatrix GetDenseDataArray(Matrix<Complex32> matrix, string name)
{
using (var stream = new MemoryStream())
using (var writer = new BinaryWriter(stream))
{
WriteMatrixTagAndName(writer, ArrayClass.Single, true, name, matrix.RowCount, matrix.ColumnCount, 0);
// write data
int size = matrix.RowCount*matrix.ColumnCount*4;
writer.Write((int)DataType.Single);
writer.Write(size);
for (var j = 0; j < matrix.ColumnCount; j++)
{ {
var column = matrix.Column(j); // small format
foreach (var value in column) writer.Write((short)DataType.Int8);
{ writer.Write((short)nameBytes.Length);
writer.Write(value.Real); writer.Write(nameBytes);
} PadData(writer, 4 - nameBytes.Length);
} }
PadData(writer, size%8); if (doublePrecision && !complex)
writer.Write((int)DataType.Single);
writer.Write(size);
for (var j = 0; j < matrix.ColumnCount; j++)
{ {
var column = matrix.Column(j); var sparseMatrix = matrix as LinearAlgebra.Double.SparseMatrix;
foreach (var value in column) if (sparseMatrix != null)
{ {
writer.Write(value.Real); SparseArrayFormatter.Write(writer, sparseMatrix);
} }
} else
PadData(writer, size%8);
writer.Flush();
return new MatlabMatrix(name, stream.ToArray());
}
}
/// <summary>
/// Gets the sparse data array.
/// </summary>
/// <param name="matrix">The matrix to get the data from.</param>
/// <param name="name">The name of the matrix.</param>
/// <returns>The matrix data as an array.</returns>
static MatlabMatrix GetSparseDataArray(LinearAlgebra.Double.SparseMatrix matrix, string name)
{
using (var stream = new MemoryStream())
using (var writer = new BinaryWriter(stream))
{
var nzmax = matrix.NonZerosCount;
WriteMatrixTagAndName(writer, ArrayClass.Sparse, false, name, matrix.RowCount, matrix.ColumnCount, nzmax);
// write ir
writer.Write((int)DataType.Int32);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{ {
writer.Write(row.Item1); NumericArrayFormatter.Write(writer, (LinearAlgebra.Double.Matrix)(object)matrix);
} }
} }
else if (!doublePrecision && !complex)
// add pad if needed
if (nzmax%2 == 1)
{ {
writer.Write(0); var sparseMatrix = matrix as LinearAlgebra.Single.SparseMatrix;
} if (sparseMatrix != null)
// write jc
writer.Write((int)DataType.Int32);
writer.Write((matrix.ColumnCount + 1)*4);
writer.Write(0);
var count = 0;
foreach (var column in matrix.EnumerateColumns())
{
count += ((SparseVectorStorage<double>)column.Storage).ValueCount;
writer.Write(count);
}
// add pad if needed
if (matrix.ColumnCount%2 == 0)
{
writer.Write(0);
}
// write data
writer.Write((int)DataType.Double);
writer.Write(nzmax*8);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{ {
writer.Write(row.Item2); SparseArrayFormatter.Write(writer, sparseMatrix);
} }
} else
writer.Flush();
return new MatlabMatrix(name, stream.ToArray());
}
}
/// <summary>
/// Gets the sparse data array.
/// </summary>
/// <param name="matrix">The matrix to get the data from.</param>
/// <param name="name">The name of the matrix.</param>
/// <returns>The matrix data as an array.</returns>
static MatlabMatrix GetSparseDataArray(LinearAlgebra.Single.SparseMatrix matrix, string name)
{
using (var stream = new MemoryStream())
using (var writer = new BinaryWriter(stream))
{
var nzmax = matrix.NonZerosCount;
WriteMatrixTagAndName(writer, ArrayClass.Sparse, false, name, matrix.RowCount, matrix.ColumnCount,
nzmax);
// write ir
writer.Write((int)DataType.Int32);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{ {
writer.Write(row.Item1); NumericArrayFormatter.Write(writer, (LinearAlgebra.Single.Matrix)(object)matrix);
} }
} }
else if (doublePrecision)
// add pad if needed
if (nzmax%2 == 1)
{ {
writer.Write(0); var sparseMatrix = matrix as LinearAlgebra.Complex.SparseMatrix;
} if (sparseMatrix != null)
// write jc
writer.Write((int)DataType.Int32);
writer.Write((matrix.ColumnCount + 1)*4);
writer.Write(0);
var count = 0;
foreach (var column in matrix.EnumerateColumns())
{
count += ((SparseVectorStorage<float>)column.Storage).ValueCount;
writer.Write(count);
}
// add pad if needed
if (matrix.ColumnCount%2 == 0)
{
writer.Write(0);
}
// write data
writer.Write((int)DataType.Single);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{ {
writer.Write(row.Item2); SparseArrayFormatter.Write(writer, sparseMatrix);
} }
} else
var pad = nzmax*4%8;
PadData(writer, pad);
writer.Flush();
return new MatlabMatrix(name, stream.ToArray());
}
}
/// <summary>
/// Gets the sparse data array.
/// </summary>
/// <param name="matrix">The matrix to get the data from.</param>
/// <param name="name">The name of the matrix.</param>
/// <returns>The matrix data as an array.</returns>
static MatlabMatrix GetSparseDataArray(LinearAlgebra.Complex.SparseMatrix matrix, string name)
{
using (var stream = new MemoryStream())
using (var writer = new BinaryWriter(stream))
{
var nzmax = matrix.NonZerosCount;
WriteMatrixTagAndName(writer, ArrayClass.Sparse, true, name, matrix.RowCount, matrix.ColumnCount,
nzmax);
// write ir
writer.Write((int)DataType.Int32);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{ {
writer.Write(row.Item1); NumericArrayFormatter.Write(writer, (LinearAlgebra.Complex.Matrix)(object)matrix);
} }
} }
else
// add pad if needed
if (nzmax%2 == 1)
{
writer.Write(0);
}
// write jc
writer.Write((int)DataType.Int32);
writer.Write((matrix.ColumnCount + 1)*4);
writer.Write(0);
var count = 0;
foreach (var column in matrix.EnumerateColumns())
{ {
count += ((SparseVectorStorage<Complex>)column.Storage).ValueCount; var sparseMatrix = matrix as LinearAlgebra.Complex32.SparseMatrix;
writer.Write(count); if (sparseMatrix != null)
}
// add pad if needed
if (matrix.ColumnCount%2 == 0)
{
writer.Write(0);
}
// write data
writer.Write((int)DataType.Double);
writer.Write(nzmax*8);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{ {
writer.Write(row.Item2.Real); SparseArrayFormatter.Write(writer, sparseMatrix);
} }
} else
writer.Write((int)DataType.Double);
writer.Write(nzmax*8);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{ {
writer.Write(row.Item2.Imaginary); NumericArrayFormatter.Write(writer, (LinearAlgebra.Complex32.Matrix)(object)matrix);
} }
} }
@ -570,105 +178,43 @@ namespace MathNet.Numerics.Data.Matlab
} }
/// <summary> /// <summary>
/// Gets the sparse data array. /// Writes all matrix blocks to a stream.
/// </summary> /// </summary>
/// <param name="matrix">The matrix to get the data from.</param> internal static void FormatFile(Stream stream, IEnumerable<MatlabMatrix> matrices)
/// <param name="name">The name of the matrix.</param>
/// <returns>The matrix data as an array.</returns>
static MatlabMatrix GetSparseDataArray(LinearAlgebra.Complex32.SparseMatrix matrix, string name)
{ {
using (var stream = new MemoryStream()) using (var buffer = new BufferedStream(stream))
using (var writer = new BinaryWriter(stream)) using (var writer = new BinaryWriter(buffer))
{ {
var nzmax = matrix.NonZerosCount; // write header and subsystem data offset (all space)
WriteMatrixTagAndName(writer, ArrayClass.Sparse, true, name, matrix.RowCount, matrix.ColumnCount, var header = Encoding.ASCII.GetBytes(HeaderText + DateTime.Now.ToString(Resources.MatlabDateHeaderFormat));
nzmax); writer.Write(header);
PadData(writer, 116 - header.Length + 8, 32);
// write ir
writer.Write((int)DataType.Int32);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{
writer.Write(row.Item1);
}
}
// add pad if needed
if (nzmax%2 == 1)
{
writer.Write(0);
}
// write jc
writer.Write((int)DataType.Int32);
writer.Write((matrix.ColumnCount + 1)*4);
writer.Write(0);
var count = 0;
foreach (var column in matrix.EnumerateColumns())
{
count += ((SparseVectorStorage<Complex32>)column.Storage).ValueCount;
writer.Write(count);
}
// add pad if needed // write version
if (matrix.ColumnCount%2 == 0) writer.Write((short)0x100);
{
writer.Write(0);
}
// write data // write little endian indicator
writer.Write((int)DataType.Single); writer.Write((byte)0x49);
writer.Write(nzmax*4); writer.Write((byte)0x4D);
foreach (var column in matrix.EnumerateColumns()) foreach (var matrix in matrices)
{ {
foreach (var row in column.EnumerateNonZeroIndexed()) // write data type
{ writer.Write((int)DataType.Compressed);
writer.Write(row.Item2.Real);
}
}
var pad = nzmax*4%8;
PadData(writer, pad);
writer.Write((int)DataType.Single); // compress data
writer.Write(nzmax*4); var compressedData = PackCompressedBlock(matrix.Data, DataType.Matrix);
foreach (var column in matrix.EnumerateColumns()) // write compressed data to file
{ writer.Write(compressedData.Length);
foreach (var row in column.EnumerateNonZeroIndexed()) writer.Write(compressedData);
{
writer.Write(row.Item2.Imaginary);
}
} }
PadData(writer, pad);
writer.Flush(); writer.Flush();
return new MatlabMatrix(name, stream.ToArray()); writer.Close();
} }
} }
/// <summary>
/// Writes the file header.
/// </summary>
static void WriteHeader(BinaryWriter writer)
{
var header = Encoding.ASCII.GetBytes(HeaderText + DateTime.Now.ToString(Resources.MatlabDateHeaderFormat));
writer.Write(header);
PadData(writer, HeaderTextLength - header.Length + 8, 32);
// write version
writer.Write((short)0x100);
// write little endian indicator
writer.Write((byte)0x49);
writer.Write((byte)0x4D);
}
/// <summary> /// <summary>
/// Pads the data with the given byte. /// Pads the data with the given byte.
/// </summary> /// </summary>
@ -684,41 +230,22 @@ namespace MathNet.Numerics.Data.Matlab
} }
/// <summary> /// <summary>
/// Writes the compressed data. /// Packs a compressed block
/// </summary> /// </summary>
/// <param name="data">The data to write.</param> static byte[] PackCompressedBlock(byte[] data, DataType dataType)
static void WriteCompressedData(BinaryWriter writer, byte[] data)
{
// fill in data size
var size = BitConverter.GetBytes(data.Length);
data[4] = size[0];
data[5] = size[1];
data[6] = size[2];
data[7] = size[3];
// compress data
var compressedData = CompressData(data);
// write compressed data to file
writer.Write(compressedData.Length);
writer.Write(compressedData);
}
/// <summary>
/// Compresses the data array.
/// </summary>
/// <param name="data">The data to compress.</param>
/// <returns>The compressed data.</returns>
static byte[] CompressData(byte[] data)
{ {
var adler = BitConverter.GetBytes(Adler32.Compute(data)); var adler = BitConverter.GetBytes(Adler32.Compute(data));
using (var compressedStream = new MemoryStream()) using (var compressedStream = new MemoryStream())
{ {
compressedStream.WriteByte(0x58); compressedStream.WriteByte(0x58);
compressedStream.WriteByte(0x85); compressedStream.WriteByte(0x85);
using (var outputStream = new DeflateStream(compressedStream, CompressionMode.Compress, true)) using (var outputStream = new DeflateStream(compressedStream, CompressionMode.Compress, true))
{ {
outputStream.Write(BitConverter.GetBytes((int)dataType), 0, 4);
outputStream.Write(BitConverter.GetBytes(data.Length), 0, 4);
outputStream.Write(data, 0, data.Length); outputStream.Write(data, 0, data.Length);
outputStream.Flush();
} }
compressedStream.WriteByte(adler[3]); compressedStream.WriteByte(adler[3]);

4
src/Data/Matlab/Matlab.csproj
View File

@ -53,13 +53,13 @@
<Compile Include="ArrayFlags.cs" /> <Compile Include="ArrayFlags.cs" />
<Compile Include="MatlabMatrix.cs" /> <Compile Include="MatlabMatrix.cs" />
<Compile Include="Formatter.cs" /> <Compile Include="Formatter.cs" />
<Compile Include="SparseArrayReader.cs" /> <Compile Include="NumericArrayFormatter.cs" />
<Compile Include="NumericArrayReader.cs" />
<Compile Include="DataType.cs" /> <Compile Include="DataType.cs" />
<Compile Include="Parser.cs" /> <Compile Include="Parser.cs" />
<Compile Include="MatlabReader.cs" /> <Compile Include="MatlabReader.cs" />
<Compile Include="MatlabWriter.cs" /> <Compile Include="MatlabWriter.cs" />
<Compile Include="Properties\AssemblyInfo.cs" /> <Compile Include="Properties\AssemblyInfo.cs" />
<Compile Include="SparseArrayFormatter.cs" />
</ItemGroup> </ItemGroup>
<ItemGroup> <ItemGroup>
<None Include="packages.config" /> <None Include="packages.config" />

6
src/Data/Matlab/MatlabMatrix.cs
View File

@ -30,14 +30,18 @@
namespace MathNet.Numerics.Data.Matlab namespace MathNet.Numerics.Data.Matlab
{ {
/// <summary>
/// MATLAB Matrix Data Element
/// </summary>
public class MatlabMatrix public class MatlabMatrix
{ {
/// <summary>Sub-elements of the matrix data element (not including the data element tag)</summary>
internal byte[] Data { get; private set; } internal byte[] Data { get; private set; }
/// <summary>Name of the matrix</summary> /// <summary>Name of the matrix</summary>
public string Name { get; private set; } public string Name { get; private set; }
/// <summary>Size of the packed matrix in bytes</summary> /// <summary>Size of the matrix in bytes</summary>
public int ByteSize public int ByteSize
{ {
get { return Data.Length; } get { return Data.Length; }

2
src/Data/Matlab/MatlabReader.cs
View File

@ -37,7 +37,7 @@ using MathNet.Numerics.LinearAlgebra;
namespace MathNet.Numerics.Data.Matlab namespace MathNet.Numerics.Data.Matlab
{ {
/// <summary> /// <summary>
/// Creates matrices from MATLAB 5 files. /// Creates matrices from MATLAB Level-5 Mat files.
/// </summary> /// </summary>
public static class MatlabReader public static class MatlabReader
{ {

2
src/Data/Matlab/MatlabWriter.cs
View File

@ -37,7 +37,7 @@ using MathNet.Numerics.LinearAlgebra;
namespace MathNet.Numerics.Data.Matlab namespace MathNet.Numerics.Data.Matlab
{ {
/// <summary> /// <summary>
/// Writes matrices to a MATLAB 5 file. /// Writes matrices to a MATLAB Level-5 Mat file.
/// </summary> /// </summary>
public static class MatlabWriter public static class MatlabWriter
{ {

150
src/Data/Matlab/NumericArrayFormatter.cs
View File

@ -0,0 +1,150 @@
// <copyright file="NumericArrayFormatter.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
// Copyright (c) 2009-2014 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
using System.IO;
using System.Numerics;
using MathNet.Numerics.LinearAlgebra;
namespace MathNet.Numerics.Data.Matlab
{
internal static class NumericArrayFormatter
{
internal static void Write(BinaryWriter writer, Matrix<double> matrix)
{
// write data
writer.Write((int)DataType.Double);
writer.Write(matrix.RowCount*matrix.ColumnCount*8);
for (var j = 0; j < matrix.ColumnCount; j++)
{
var column = matrix.Column(j);
foreach (var value in column)
{
writer.Write(value);
}
}
}
internal static void Write(BinaryWriter writer, Matrix<float> matrix)
{
// write data
int size = matrix.RowCount*matrix.ColumnCount*4;
writer.Write((int)DataType.Single);
writer.Write(size);
for (var j = 0; j < matrix.ColumnCount; j++)
{
var column = matrix.Column(j);
foreach (var value in column)
{
writer.Write(value);
}
}
PadData(writer, size%8);
}
internal static void Write(BinaryWriter writer, Matrix<Complex> matrix)
{
// write data
int size = matrix.RowCount*matrix.ColumnCount*8;
writer.Write((int)DataType.Double);
writer.Write(size);
for (var j = 0; j < matrix.ColumnCount; j++)
{
var column = matrix.Column(j);
foreach (var value in column)
{
writer.Write(value.Real);
}
}
writer.Write((int)DataType.Double);
writer.Write(size);
for (var j = 0; j < matrix.ColumnCount; j++)
{
var column = matrix.Column(j);
foreach (var value in column)
{
writer.Write(value.Imaginary);
}
}
}
internal static void Write(BinaryWriter writer, Matrix<Complex32> matrix)
{
// write data
int size = matrix.RowCount*matrix.ColumnCount*4;
writer.Write((int)DataType.Single);
writer.Write(size);
for (var j = 0; j < matrix.ColumnCount; j++)
{
var column = matrix.Column(j);
foreach (var value in column)
{
writer.Write(value.Real);
}
}
PadData(writer, size%8);
writer.Write((int)DataType.Single);
writer.Write(size);
for (var j = 0; j < matrix.ColumnCount; j++)
{
var column = matrix.Column(j);
foreach (var value in column)
{
writer.Write(value.Real);
}
}
PadData(writer, size%8);
}
/// <summary>
/// Pads the data with the given byte.
/// </summary>
/// <param name="writer">Where to write the pad values.</param>
/// <param name="bytes">The number of bytes to pad.</param>
/// <param name="pad">What value to pad with.</param>
static void PadData(BinaryWriter writer, int bytes, byte pad = (byte)0)
{
for (var i = 0; i < bytes; i++)
{
writer.Write(pad);
}
}
}
}

257
src/Data/Matlab/NumericArrayReader.cs
View File

@ -1,257 +0,0 @@
// <copyright file="NumericArrayReader.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
// Copyright (c) 2009-2014 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
using System;
using System.IO;
using System.Numerics;
using MathNet.Numerics.LinearAlgebra;
namespace MathNet.Numerics.Data.Matlab
{
internal static class NumericArrayReader<TDataType>
where TDataType : struct, IEquatable<TDataType>, IFormattable
{
/// <summary>
/// Populates a dense matrix.
/// </summary>
/// <param name="type">The MATLAB data type.</param>
/// <param name="reader">The reader to read from.</param>
/// <param name="isComplex">if set to <c>true</c> if the MATLAB complex flag is set.</param>
/// <param name="rows">The number of rows.</param>
/// <param name="columns">The number of columns.</param>
/// <param name="size">The length of the stored data.</param>
/// <returns>Returns a populated dense matrix.</returns>
public static Matrix<TDataType> PopulateDenseMatrix(DataType type, BinaryReader reader, bool isComplex, int rows, int columns, int size)
{
var dataType = typeof (TDataType);
Matrix<TDataType> matrix;
if (type == DataType.Double && dataType == typeof (double))
{
var count = rows*columns;
var data = new double[count];
Buffer.BlockCopy(reader.ReadBytes(count*Constants.SizeOfDouble), 0, data, 0, count*Constants.SizeOfDouble);
matrix = (Matrix<TDataType>)(object)new LinearAlgebra.Double.DenseMatrix(rows, columns, data);
}
else if (type == DataType.Single && dataType == typeof (float))
{
var count = rows*columns;
var data = new float[count];
Buffer.BlockCopy(reader.ReadBytes(count*Constants.SizeOfFloat), 0, data, 0, count*Constants.SizeOfFloat);
matrix = (Matrix<TDataType>)(object)new LinearAlgebra.Single.DenseMatrix(rows, columns, data);
}
else
{
matrix = Matrix<TDataType>.Build.Dense(rows, columns);
if (dataType == typeof (double))
{
if (isComplex)
{
throw new ArgumentException("Invalid TDataType. Matrix is stored as a complex matrix, but a real data type was given.");
}
PopulateDoubleDenseMatrix((Matrix<double>)(object)matrix, type, reader, rows, columns);
}
else if (dataType == typeof (float))
{
if (isComplex)
{
throw new ArgumentException("Invalid TDataType. Matrix is stored as a complex matrix, but a real data type was given.");
}
PopulateSingleDenseMatrix((Matrix<float>)(object)matrix, type, reader, rows, columns);
}
else if (dataType == typeof (Complex))
{
PopulateComplexDenseMatrix((Matrix<Complex>)(object)matrix, type, isComplex, reader, rows, columns, size);
}
else if (dataType == typeof (Complex32))
{
PopulateComplex32DenseMatrix((Matrix<Complex32>)(object)matrix, type, isComplex, reader, rows, columns, size);
}
else
{
throw new NotSupportedException();
}
}
return matrix;
}
/// <summary>
/// Populates the double dense matrix.
/// </summary>
/// <param name="matrix">The matrix to populate.</param>
/// <param name="type">The MATLAB data type.</param>
/// <param name="reader">The reader to read from.</param>
/// <param name="rows">The number of rows.</param>
/// <param name="columns">The number of columns.</param>
static void PopulateDoubleDenseMatrix(Matrix<double> matrix, DataType type, BinaryReader reader, int rows, int columns)
{
for (var j = 0; j < columns; j++)
{
for (var i = 0; i < rows; i++)
{
matrix.At(i, j, ReadDoubleValue(type, reader));
}
}
}
/// <summary>
/// Populates the complex dense matrix.
/// </summary>
/// <param name="matrix">The matrix to populate.</param>
/// <param name="type">The MATLAB data type.</param>
/// <param name="isComplex">if set to <c>true</c> if the MATLAB complex flag is set.</param>
/// <param name="reader">The reader to read from.</param>
/// <param name="rows">The number of rows.</param>
/// <param name="columns">The number of columns.</param>
/// <param name="dataSize">The length of the stored data.</param>
static void PopulateComplexDenseMatrix(Matrix<Complex> matrix, DataType type, bool isComplex, BinaryReader reader, int rows, int columns, int dataSize)
{
for (var j = 0; j < columns; j++)
{
for (var i = 0; i < rows; i++)
{
matrix.At(i, j, ReadDoubleValue(type, reader));
}
}
if (isComplex)
{
var skip = dataSize%8;
// skip pad
reader.ReadBytes(skip);
// skip header
type = (DataType)reader.ReadInt32();
reader.ReadInt32();
for (var j = 0; j < columns; j++)
{
for (var i = 0; i < rows; i++)
{
matrix.At(i, j, new Complex(matrix.At(i, j).Real, ReadDoubleValue(type, reader)));
}
}
}
}
/// <summary>
/// Populates the complex32 dense matrix.
/// </summary>
/// <param name="matrix">The matrix to populate.</param>
/// <param name="type">The MATLAB data type.</param>
/// <param name="isComplex">if set to <c>true</c> if the MATLAB complex flag is set.</param>
/// <param name="reader">The reader to read from.</param>
/// <param name="rows">The number of rows.</param>
/// <param name="columns">The number of columns.</param>
/// <param name="dataSize">The length of the stored data.</param>
static void PopulateComplex32DenseMatrix(Matrix<Complex32> matrix, DataType type, bool isComplex, BinaryReader reader, int rows, int columns, int dataSize)
{
for (var j = 0; j < columns; j++)
{
for (var i = 0; i < rows; i++)
{
matrix.At(i, j, (float)ReadDoubleValue(type, reader));
}
}
if (isComplex)
{
var skip = dataSize%8;
// skip pad
reader.ReadBytes(skip);
// skip header
type = (DataType)reader.ReadInt32();
reader.ReadInt32();
for (var j = 0; j < columns; j++)
{
for (var i = 0; i < rows; i++)
{
matrix.At(i, j, new Complex32(matrix.At(i, j).Real, (float)ReadDoubleValue(type, reader)));
}
}
}
}
/// <summary>
/// Populates the float dense matrix.
/// </summary>
/// <param name="matrix">The matrix to populate.</param>
/// <param name="type">The MATLAB data type.</param>
/// <param name="reader">The reader to read from.</param>
/// <param name="rows">The number of rows.</param>
/// <param name="columns">The number of columns.</param>
static void PopulateSingleDenseMatrix(Matrix<float> matrix, DataType type, BinaryReader reader, int rows, int columns)
{
for (var j = 0; j < columns; j++)
{
for (var i = 0; i < rows; i++)
{
matrix.At(i, j, (float)ReadDoubleValue(type, reader));
}
}
}
static double ReadDoubleValue(DataType type, BinaryReader reader)
{
switch (type)
{
case DataType.Double:
return reader.ReadDouble();
case DataType.Int8:
return reader.ReadSByte();
case DataType.UInt8:
return reader.ReadByte();
case DataType.Int16:
return reader.ReadInt16();
case DataType.UInt16:
return reader.ReadUInt16();
case DataType.Int32:
return reader.ReadInt32();
case DataType.UInt32:
return reader.ReadUInt32();
case DataType.Single:
return reader.ReadSingle();
case DataType.Int64:
return reader.ReadInt64();
case DataType.UInt64:
return reader.ReadUInt64();
default:
throw new NotSupportedException();
}
}
}
}

447
src/Data/Matlab/Parser.cs
View File

@ -32,8 +32,10 @@ using System;
using System.Collections.Generic; using System.Collections.Generic;
using System.IO; using System.IO;
using System.IO.Compression; using System.IO.Compression;
using System.Numerics;
using System.Text; using System.Text;
using MathNet.Numerics.LinearAlgebra; using MathNet.Numerics.LinearAlgebra;
using MathNet.Numerics.LinearAlgebra.Storage;
using MathNet.Numerics.Properties; using MathNet.Numerics.Properties;
namespace MathNet.Numerics.Data.Matlab namespace MathNet.Numerics.Data.Matlab
@ -43,11 +45,6 @@ namespace MathNet.Numerics.Data.Matlab
/// </summary> /// </summary>
internal static class Parser internal static class Parser
{ {
/// <summary>
/// Large Block Size
/// </summary>
const int LargeBlockSize = 8;
/// <summary> /// <summary>
/// Little Endian Indicator /// Little Endian Indicator
/// </summary> /// </summary>
@ -59,74 +56,9 @@ namespace MathNet.Numerics.Data.Matlab
const int SmallBlockSize = 4; const int SmallBlockSize = 4;
/// <summary> /// <summary>
/// Parse a matrix block byte array /// Large Block Size
/// </summary> /// </summary>
internal static Matrix<T> ParseMatrix<T>(byte[] data) const int LargeBlockSize = 8;
where T : struct, IEquatable<T>, IFormattable
{
using (var stream = new MemoryStream(data))
using (var reader = new BinaryReader(stream))
{
// skip tag - doesn't tell us anything we don't already know
reader.BaseStream.Seek(8, SeekOrigin.Current);
var arrayClass = (ArrayClass)reader.ReadByte();
var flags = reader.ReadByte();
var isComplex = (flags & (byte)ArrayFlags.Complex) == (byte)ArrayFlags.Complex;
// skip unneeded bytes
reader.BaseStream.Seek(10, SeekOrigin.Current);
var numDimensions = reader.ReadInt32()/8;
if (numDimensions > 2)
{
throw new NotSupportedException(Resources.MoreThan2D);
}
var rows = reader.ReadInt32();
var columns = reader.ReadInt32();
// skip name and unneeded bytes
reader.BaseStream.Seek(2, SeekOrigin.Current);
int size = reader.ReadInt16();
var smallBlock = true;
if (size == 0)
{
size = reader.ReadInt32();
smallBlock = false;
}
reader.BaseStream.Seek(size, SeekOrigin.Current);
AlignData(reader.BaseStream, size, smallBlock);
var type = (DataType)reader.ReadInt16();
size = reader.ReadInt16();
if (size == 0)
{
size = reader.ReadInt32();
}
Matrix<T> matrix;
switch (arrayClass)
{
case ArrayClass.Sparse:
matrix = SparseArrayReader<T>.PopulateSparseMatrix(reader, isComplex, rows, columns, size);
break;
case ArrayClass.Function:
case ArrayClass.Character:
case ArrayClass.Object:
case ArrayClass.Structure:
case ArrayClass.Cell:
case ArrayClass.Unknown:
throw new NotSupportedException();
default:
matrix = NumericArrayReader<T>.PopulateDenseMatrix(type, reader, isComplex, rows, columns, size);
break;
}
return matrix;
}
}
/// <summary> /// <summary>
/// Extracts all matrix blocks in a format we support from a stream. /// Extracts all matrix blocks in a format we support from a stream.
@ -137,38 +69,44 @@ namespace MathNet.Numerics.Data.Matlab
using (var reader = new BinaryReader(stream)) using (var reader = new BinaryReader(stream))
{ {
// skip header (116 bytes)
// skip subsystem data offset (8 bytes)
// skip version (2 bytes)
reader.BaseStream.Position = 126; reader.BaseStream.Position = 126;
// endian indicator (2 bytes)
if (reader.ReadByte() != LittleEndianIndicator) if (reader.ReadByte() != LittleEndianIndicator)
{ {
throw new NotSupportedException(Resources.BigEndianNotSupported); throw new NotSupportedException(Resources.BigEndianNotSupported);
} }
// skip version since it is always 0x0100. // set position to first data element, right after full file header (128 bytes)
reader.BaseStream.Position = 128; reader.BaseStream.Position = 128;
var length = stream.Length; var length = stream.Length;
// for each data block add a MATLAB object to the file. // for each data element add a MATLAB object to the file.
while (reader.BaseStream.Position < length) while (reader.BaseStream.Position < length)
{ {
var type = (DataType)reader.ReadInt16(); // small format: size (2 bytes), type (2 bytes), data (4 bytes)
int size = reader.ReadInt16(); // long format: type (4 bytes), size (4 bytes), data (size, aligned to 8 bytes)
var smallBlock = true;
if (size == 0)
{
size = reader.ReadInt32();
smallBlock = false;
}
DataType type;
int size;
bool smallBlock;
ReadElementTag(reader, out type, out size, out smallBlock);
// read element data of the size provided in the element header
// uncompress if compressed
byte[] data; byte[] data;
if (type == DataType.Compressed) if (type == DataType.Compressed)
{ {
data = DecompressBlock(reader.ReadBytes(size), out type); data = UnpackCompressedBlock(reader.ReadBytes(size), out type);
} }
else else
{ {
data = new byte[size]; data = new byte[size];
reader.Read(data, 0, size); reader.Read(data, 0, size);
AlignData(reader.BaseStream, size, smallBlock); SkipElementPadding(reader.BaseStream, size, smallBlock);
} }
if (type == DataType.Matrix) if (type == DataType.Matrix)
@ -202,12 +140,330 @@ namespace MathNet.Numerics.Data.Matlab
} }
/// <summary> /// <summary>
/// Aligns the data. /// Parse a matrix block byte array
/// </summary>
internal static Matrix<T> ParseMatrix<T>(byte[] data)
where T : struct, IEquatable<T>, IFormattable
{
using (var stream = new MemoryStream(data))
using (var reader = new BinaryReader(stream))
{
// Array Flags tag (8 bytes)
reader.BaseStream.Seek(8, SeekOrigin.Current);
// Array Flags data: flags (byte 3), class (byte 4) (8 bytes)
var arrayClass = (ArrayClass)reader.ReadByte();
var flags = reader.ReadByte();
var complex = (flags & (byte)ArrayFlags.Complex) == (byte)ArrayFlags.Complex;
reader.BaseStream.Seek(6, SeekOrigin.Current);
// Dimensions Array tag (8 bytes)
reader.BaseStream.Seek(4, SeekOrigin.Current);
var numDimensions = reader.ReadInt32()/8;
if (numDimensions > 2)
{
throw new NotSupportedException(Resources.MoreThan2D);
}
// Dimensions Array data: row and column count (8 bytes)
var rows = reader.ReadInt32();
var columns = reader.ReadInt32();
// Array name
DataType type;
int size;
bool smallBlock;
ReadElementTag(reader, out type, out size, out smallBlock);
reader.BaseStream.Seek(size, SeekOrigin.Current);
SkipElementPadding(reader.BaseStream, size, smallBlock);
Matrix<T> matrix;
switch (arrayClass)
{
case ArrayClass.Sparse:
matrix = PopulateSparseMatrix<T>(reader, complex, rows, columns);
break;
case ArrayClass.Function:
case ArrayClass.Character:
case ArrayClass.Object:
case ArrayClass.Structure:
case ArrayClass.Cell:
case ArrayClass.Unknown:
throw new NotSupportedException();
default:
matrix = PopulateDenseMatrix<T>(reader, complex, rows, columns);
break;
}
return matrix;
}
}
/// <summary>
/// Populates a dense matrix.
/// </summary>
/// <param name="reader">The reader to read from.</param>
/// <param name="complex">if set to <c>true</c> if the MATLAB complex flag is set.</param>
/// <param name="rows">The number of rows.</param>
/// <param name="columns">The number of columns.</param>
/// <returns>Returns a populated dense matrix.</returns>
static Matrix<T> PopulateDenseMatrix<T>(BinaryReader reader, bool complex, int rows, int columns)
where T : struct, IEquatable<T>, IFormattable
{
var dataType = typeof(T);
var count = rows*columns;
var data = new T[count];
DataType type;
int size;
bool smallBlock;
// read real part array
ReadElementTag(reader, out type, out size, out smallBlock);
// direct copy if possible
if (type == DataType.Double && dataType == typeof(double) || type == DataType.Single && dataType == typeof(float))
{
Buffer.BlockCopy(reader.ReadBytes(size), 0, data, 0, size);
}
else if (dataType == typeof(double))
{
if (complex)
{
throw new ArgumentException("Invalid TDataType. Matrix is stored as a complex matrix, but a real data type was given.");
}
PopulateDoubleArray(reader, (double[])(object)data, type);
}
else if (dataType == typeof(float))
{
if (complex)
{
throw new ArgumentException("Invalid TDataType. Matrix is stored as a complex matrix, but a real data type was given.");
}
PopulateSingleArray(reader, (float[])(object)data, type);
}
else if (dataType == typeof(Complex))
{
PopulateComplexArray(reader, (Complex[])(object)data, complex, type, ref size, ref smallBlock);
}
else if (dataType == typeof(Complex32))
{
PopulateComplex32Array(reader, (Complex32[])(object)data, complex, type, ref size, ref smallBlock);
}
else
{
throw new NotSupportedException();
}
SkipElementPadding(reader.BaseStream, size, smallBlock);
return Matrix<T>.Build.Dense(rows, columns, data);
}
/// <summary>
/// Populates a sparse matrix.
/// </summary> /// </summary>
/// <param name="stream">The stream.</param> /// <param name="reader">The reader.</param>
/// <param name="size">The size of the array.</param> /// <param name="complex">if set to <c>true</c> if the MATLAB complex flag is set.</param>
/// <param name="smallBlock">if set to <c>true</c> if reading from a small block.</param> /// <param name="rows">The number of rows.</param>
internal static void AlignData(Stream stream, int size, bool smallBlock) /// <param name="columns">The number of columns.</param>
/// <returns>A populated sparse matrix.</returns>
static Matrix<T> PopulateSparseMatrix<T>(BinaryReader reader, bool complex, int rows, int columns)
where T : struct, IEquatable<T>, IFormattable
{
// Create matrix with CSR storage.
var matrix = Matrix<T>.Build.Sparse(columns, rows);
// MATLAB sparse matrices are actually stored as CSC, so just read the data and then transpose.
var storage = matrix.Storage as SparseCompressedRowMatrixStorage<T>;
DataType type;
int size;
bool smallBlock;
// populate the row data array
ReadElementTag(reader, out type, out size, out smallBlock);
var ir = storage.ColumnIndices = new int[size/4];
for (var i = 0; i < ir.Length; i++)
{
ir[i] = reader.ReadInt32();
}
SkipElementPadding(reader.BaseStream, size, smallBlock);
// populate the column data array
ReadElementTag(reader, out type, out size, out smallBlock);
var jc = storage.RowPointers;
if (jc.Length != size/4)
{
throw new Exception("invalid jcsize");
}
for (var j = 0; j < jc.Length; j++)
{
jc[j] = reader.ReadInt32();
}
SkipElementPadding(reader.BaseStream, size, smallBlock);
// populate the values
ReadElementTag(reader, out type, out size, out smallBlock);
var dataType = typeof(T);
var data = storage.Values = new T[jc[columns]];
if (dataType == typeof(double))
{
if (complex)
{
throw new ArgumentException("Invalid TDataType. Matrix is stored as a complex matrix, but a real data type was given.");
}
PopulateDoubleArray(reader, (double[])(object)data, type);
}
else if (dataType == typeof(float))
{
if (complex)
{
throw new ArgumentException("Invalid TDataType. Matrix is stored as a complex matrix, but a real data type was given.");
}
PopulateSingleArray(reader, (float[])(object)data, type);
}
else if (dataType == typeof(Complex))
{
PopulateComplexArray(reader, (Complex[])(object)data, complex, type, ref size, ref smallBlock);
}
else if (dataType == typeof(Complex32))
{
PopulateComplex32Array(reader, (Complex32[])(object)data, complex, type, ref size, ref smallBlock);
}
else
{
throw new NotSupportedException();
}
SkipElementPadding(reader.BaseStream, size, smallBlock);
return matrix.Transpose();
}
/// <summary>
/// Populates the double dense matrix.
/// </summary>
static void PopulateDoubleArray(BinaryReader reader, double[] data, DataType type)
{
for (int i = 0; i < data.Length; i++)
{
data[i] = ReadDoubleValue(reader, type);
}
}
/// <summary>
/// Populates the float dense matrix.
/// </summary>
static void PopulateSingleArray(BinaryReader reader, float[] data, DataType type)
{
for (int i = 0; i < data.Length; i++)
{
data[i] = (float)ReadDoubleValue(reader, type);
}
}
/// <summary>
/// Populates the complex dense matrix.
/// </summary>
static void PopulateComplexArray(BinaryReader reader, Complex[] data, bool complex, DataType type, ref int size, ref bool smallBlock)
{
for (int i = 0; i < data.Length; i++)
{
data[i] = ReadDoubleValue(reader, type);
}
if (complex)
{
SkipElementPadding(reader.BaseStream, size, smallBlock);
ReadElementTag(reader, out type, out size, out smallBlock);
for (int i = 0; i < data.Length; i++)
{
data[i] = new Complex(data[i].Real, ReadDoubleValue(reader, type));
}
}
}
/// <summary>
/// Populates the complex32 dense matrix.
/// </summary>
static void PopulateComplex32Array(BinaryReader reader, Complex32[] data, bool complex, DataType type, ref int size, ref bool smallBlock)
{
for (int i = 0; i < data.Length; i++)
{
data[i] = (float)ReadDoubleValue(reader, type);
}
if (complex)
{
SkipElementPadding(reader.BaseStream, size, smallBlock);
ReadElementTag(reader, out type, out size, out smallBlock);
for (int i = 0; i < data.Length; i++)
{
data[i] = new Complex32(data[i].Real, (float)ReadDoubleValue(reader, type));
}
}
}
static double ReadDoubleValue(BinaryReader reader, DataType type)
{
switch (type)
{
case DataType.Double:
return reader.ReadDouble();
case DataType.Int8:
return reader.ReadSByte();
case DataType.UInt8:
return reader.ReadByte();
case DataType.Int16:
return reader.ReadInt16();
case DataType.UInt16:
return reader.ReadUInt16();
case DataType.Int32:
return reader.ReadInt32();
case DataType.UInt32:
return reader.ReadUInt32();
case DataType.Single:
return reader.ReadSingle();
case DataType.Int64:
return reader.ReadInt64();
case DataType.UInt64:
return reader.ReadUInt64();
default:
throw new NotSupportedException();
}
}
static void ReadElementTag(BinaryReader reader, out DataType dataType, out int size, out bool smallBlock)
{
// assume small format
smallBlock = true;
// small type (2 bytes)
dataType = (DataType)reader.ReadInt16();
// small size (2 bytes)
size = reader.ReadInt16();
if (size == 0)
{
// long format detected
smallBlock = false;
// long size (4 bytes)
size = reader.ReadInt32();
}
}
static void SkipElementPadding(Stream stream, int size, bool smallBlock)
{ {
var blockSize = smallBlock ? SmallBlockSize : LargeBlockSize; var blockSize = smallBlock ? SmallBlockSize : LargeBlockSize;
var offset = 0; var offset = 0;
@ -221,19 +477,22 @@ namespace MathNet.Numerics.Data.Matlab
} }
/// <summary> /// <summary>
/// Decompresses the block. /// Unpacks a compressed block.
/// </summary> /// </summary>
/// <param name="compressed">The compressed data.</param> /// <param name="compressed">The compressed data.</param>
/// <param name="type">The type data type contained in the block.</param> /// <param name="type">The type data type contained in the block.</param>
/// <returns>The decompressed block.</returns> /// <returns>The decompressed block.</returns>
static byte[] DecompressBlock(byte[] compressed, out DataType type) static byte[] UnpackCompressedBlock(byte[] compressed, out DataType type)
{ {
byte[] data; byte[] data;
using (var compressedStream = new MemoryStream(compressed, 2, compressed.Length - 6))
using (var decompressor = new DeflateStream(compressedStream, CompressionMode.Decompress))
using (var decompressed = new MemoryStream()) using (var decompressed = new MemoryStream())
{ {
decompressor.CopyTo(decompressed); using (var compressedStream = new MemoryStream(compressed, 2, compressed.Length - 6))
using (var decompressor = new DeflateStream(compressedStream, CompressionMode.Decompress))
{
decompressor.CopyTo(decompressed);
}
decompressed.Position = 0; decompressed.Position = 0;
var buf = new byte[4]; var buf = new byte[4];
decompressed.Read(buf, 0, 4); decompressed.Read(buf, 0, 4);

301
src/Data/Matlab/SparseArrayFormatter.cs
View File

@ -0,0 +1,301 @@
// <copyright file="SparseArrayFormatter.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
// Copyright (c) 2009-2014 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
using System.IO;
using System.Numerics;
using MathNet.Numerics.LinearAlgebra.Storage;
namespace MathNet.Numerics.Data.Matlab
{
internal static class SparseArrayFormatter
{
//public static void FormatSparseMatrix<T>(BinaryWriter writer, Matrix<T> matrix, string name, bool isComplex, int rows, int columns, int size)
// where T : struct, IEquatable<T>, IFormattable
//{
// var transposed = matrix.Transpose();
// var storage = (SparseCompressedRowMatrixStorage<T>)transposed.Storage;
// WriteMatrixTagAndName(writer, ArrayClass.Sparse, isComplex, name, storage);
//}
internal static void Write(BinaryWriter writer, LinearAlgebra.Double.SparseMatrix matrix)
{
var nzmax = matrix.NonZerosCount;
// write ir
writer.Write((int)DataType.Int32);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{
writer.Write(row.Item1);
}
}
// add pad if needed
if (nzmax%2 == 1)
{
writer.Write(0);
}
// write jc
writer.Write((int)DataType.Int32);
writer.Write((matrix.ColumnCount + 1)*4);
writer.Write(0);
var count = 0;
foreach (var column in matrix.EnumerateColumns())
{
count += ((SparseVectorStorage<double>)column.Storage).ValueCount;
writer.Write(count);
}
// add pad if needed
if (matrix.ColumnCount%2 == 0)
{
writer.Write(0);
}
// write data
writer.Write((int)DataType.Double);
writer.Write(nzmax*8);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{
writer.Write(row.Item2);
}
}
}
internal static void Write(BinaryWriter writer, LinearAlgebra.Single.SparseMatrix matrix)
{
var nzmax = matrix.NonZerosCount;
// write ir
writer.Write((int)DataType.Int32);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{
writer.Write(row.Item1);
}
}
// add pad if needed
if (nzmax%2 == 1)
{
writer.Write(0);
}
// write jc
writer.Write((int)DataType.Int32);
writer.Write((matrix.ColumnCount + 1)*4);
writer.Write(0);
var count = 0;
foreach (var column in matrix.EnumerateColumns())
{
count += ((SparseVectorStorage<float>)column.Storage).ValueCount;
writer.Write(count);
}
// add pad if needed
if (matrix.ColumnCount%2 == 0)
{
writer.Write(0);
}
// write data
writer.Write((int)DataType.Single);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{
writer.Write(row.Item2);
}
}
var pad = nzmax*4%8;
PadData(writer, pad);
}
internal static void Write(BinaryWriter writer, LinearAlgebra.Complex.SparseMatrix matrix)
{
var nzmax = matrix.NonZerosCount;
// write ir
writer.Write((int)DataType.Int32);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{
writer.Write(row.Item1);
}
}
// add pad if needed
if (nzmax%2 == 1)
{
writer.Write(0);
}
// write jc
writer.Write((int)DataType.Int32);
writer.Write((matrix.ColumnCount + 1)*4);
writer.Write(0);
var count = 0;
foreach (var column in matrix.EnumerateColumns())
{
count += ((SparseVectorStorage<Complex>)column.Storage).ValueCount;
writer.Write(count);
}
// add pad if needed
if (matrix.ColumnCount%2 == 0)
{
writer.Write(0);
}
// write data
writer.Write((int)DataType.Double);
writer.Write(nzmax*8);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{
writer.Write(row.Item2.Real);
}
}
writer.Write((int)DataType.Double);
writer.Write(nzmax*8);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{
writer.Write(row.Item2.Imaginary);
}
}
}
internal static void Write(BinaryWriter writer, LinearAlgebra.Complex32.SparseMatrix matrix)
{
var nzmax = matrix.NonZerosCount;
// write ir
writer.Write((int)DataType.Int32);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{
writer.Write(row.Item1);
}
}
// add pad if needed
if (nzmax%2 == 1)
{
writer.Write(0);
}
// write jc
writer.Write((int)DataType.Int32);
writer.Write((matrix.ColumnCount + 1)*4);
writer.Write(0);
var count = 0;
foreach (var column in matrix.EnumerateColumns())
{
count += ((SparseVectorStorage<Complex32>)column.Storage).ValueCount;
writer.Write(count);
}
// add pad if needed
if (matrix.ColumnCount%2 == 0)
{
writer.Write(0);
}
// write data
writer.Write((int)DataType.Single);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{
writer.Write(row.Item2.Real);
}
}
var pad = nzmax*4%8;
PadData(writer, pad);
writer.Write((int)DataType.Single);
writer.Write(nzmax*4);
foreach (var column in matrix.EnumerateColumns())
{
foreach (var row in column.EnumerateNonZeroIndexed())
{
writer.Write(row.Item2.Imaginary);
}
}
PadData(writer, pad);
}
/// <summary>
/// Pads the data with the given byte.
/// </summary>
/// <param name="writer">Where to write the pad values.</param>
/// <param name="bytes">The number of bytes to pad.</param>
/// <param name="pad">What value to pad with.</param>
static void PadData(BinaryWriter writer, int bytes, byte pad = (byte)0)
{
for (var i = 0; i < bytes; i++)
{
writer.Write(pad);
}
}
}
}

251
src/Data/Matlab/SparseArrayReader.cs
View File

@ -1,251 +0,0 @@
// <copyright file="SparseArrayReader.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
// http://mathnetnumerics.codeplex.com
//
// Copyright (c) 2009-2014 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
using System;
using System.IO;
using System.Numerics;
using MathNet.Numerics.LinearAlgebra;
using MathNet.Numerics.LinearAlgebra.Storage;
namespace MathNet.Numerics.Data.Matlab
{
internal static class SparseArrayReader<TDataType>
where TDataType : struct, IEquatable<TDataType>, IFormattable
{
/// <summary>
/// Populates a sparse matrix.
/// </summary>
/// <param name="reader">The reader.</param>
/// <param name="isComplex">if set to <c>true</c> if the MATLAB complex flag is set.</param>
/// <param name="rows">The number of rows.</param>
/// <param name="columns">The number of columns.</param>
/// <param name="size">The size of the block.</param>
/// <returns>A populated sparse matrix.</returns>
public static Matrix<TDataType> PopulateSparseMatrix(BinaryReader reader, bool isComplex, int rows, int columns, int size)
{
// Create matrix with CSR storage.
var matrix = Matrix<TDataType>.Build.Sparse(columns, rows);
// MATLAB sparse matrices are actually stored as CSC, so just read the data and then transpose.
var storage = matrix.Storage as SparseCompressedRowMatrixStorage<TDataType>;
// populate the row data array
var ir = storage.ColumnIndices = new int[size/4];
for (var i = 0; i < ir.Length; i++)
{
ir[i] = reader.ReadInt32();
}
Parser.AlignData(reader.BaseStream, size, false);
// skip data type since it will always be int32
reader.BaseStream.Seek(4, SeekOrigin.Current);
// populate the column data array
var jcsize = reader.ReadInt32();
var jc = storage.RowPointers;
if (jc.Length != jcsize/4)
{
throw new Exception("invalid jcsize");
}
for (var j = 0; j < jc.Length; j++)
{
jc[j] = reader.ReadInt32();
}
Parser.AlignData(reader.BaseStream, jcsize, false);
var type = (DataType)reader.ReadInt32();
var dataSize = reader.ReadInt32();
var dataType = typeof (TDataType);
// Allocate memory for matrix values
var data = storage.Values = new TDataType[jc[columns]];
if (dataType == typeof (double))
{
if (isComplex)
{
throw new ArgumentException("Invalid TDataType. Matrix is stored as a complex matrix, but a real data type was given.");
}
PopulateDoubleSparseMatrix(type, (double[])(object)data, reader);
}
else if (dataType == typeof (float))
{
if (isComplex)
{
throw new ArgumentException("Invalid TDataType. Matrix is stored as a complex matrix, but a real data type was given.");
}
PopulateSingleSparseMatrix(type, (float[])(object)data, reader);
}
else if (dataType == typeof (Complex))
{
PopulateComplexSparseMatrix(type, isComplex, (Complex[])(object)data, reader, dataSize);
}
else if (dataType == typeof (Complex32))
{
PopulateComplex32SparseMatrix(type, isComplex, (Complex32[])(object)data, reader, dataSize);
}
else
{
throw new NotSupportedException();
}
return matrix.Transpose();
}
/// <summary>
/// Populates the double sparse matrix.
/// </summary>
/// <param name="type">The MATLAB data type.</param>
/// <param name="data">The matrix values array.</param>
/// <param name="reader">The reader to read from.</param>
static void PopulateDoubleSparseMatrix(DataType type, double[] data, BinaryReader reader)
{
for (var i = 0; i < data.Length; i++)
{
data[i] = ReadDoubleValue(type, reader);
}
}
/// <summary>
/// Populates the float sparse matrix.
/// </summary>
/// <param name="type">The MATLAB data type.</param>
/// <param name="data">The matrix values array.</param>
/// <param name="reader">The reader to read from.</param>
static void PopulateSingleSparseMatrix(DataType type, float[] data, BinaryReader reader)
{
for (var i = 0; i < data.Length; i++)
{
data[i] = (float)ReadDoubleValue(type, reader);
}
}
/// <summary>
/// Populates the complex sparse matrix.
/// </summary>
/// <param name="type">The MATLAB data type.</param>
/// <param name="isComplex">if set to <c>true</c> if the MATLAB complex flag is set.</param>
/// <param name="data">The matrix values array.</param>
/// <param name="reader">The reader to read from.</param>
/// <param name="dataSize">The length of the stored data.</param>
static void PopulateComplexSparseMatrix(DataType type, bool isComplex, Complex[] data, BinaryReader reader, int dataSize)
{
for (var i = 0; i < data.Length; i++)
{
data[i] = ReadDoubleValue(type, reader);
}
if (isComplex)
{
var skip = dataSize%8;
// skip pad
reader.ReadBytes(skip);
// skip header
type = (DataType)reader.ReadInt32();
reader.ReadInt32();
for (var i = 0; i < data.Length; i++)
{
data[i] += new Complex(0.0, ReadDoubleValue(type, reader));
}
}
}
/// <summary>
/// Populates the complex32 sparse matrix.
/// </summary>
/// <param name="type">The MATLAB data type.</param>
/// <param name="isComplex">if set to <c>true</c> if the MATLAB complex flag is set.</param>
/// <param name="data">The matrix values array.</param>
/// <param name="reader">The reader to read from.</param>
/// <param name="dataSize">The length of the stored data.</param>
static void PopulateComplex32SparseMatrix(DataType type, bool isComplex, Complex32[] data, BinaryReader reader, int dataSize)
{
for (var i = 0; i < data.Length; i++)
{
data[i] = (float)ReadDoubleValue(type, reader);
}
if (isComplex)
{
var skip = dataSize%8;
// skip pad
reader.ReadBytes(skip);
// skip header
type = (DataType)reader.ReadInt32();
reader.ReadInt32();
for (var i = 0; i < data.Length; i++)
{
data[i] += new Complex32(0.0f, (float)ReadDoubleValue(type, reader));
}
}
}
static double ReadDoubleValue(DataType type, BinaryReader reader)
{
switch (type)
{
case DataType.Double:
return reader.ReadDouble();
case DataType.Int8:
return reader.ReadSByte();
case DataType.UInt8:
return reader.ReadByte();
case DataType.Int16:
return reader.ReadInt16();
case DataType.UInt16:
return reader.ReadUInt16();
case DataType.Int32:
return reader.ReadInt32();
case DataType.UInt32:
return reader.ReadUInt32();
case DataType.Single:
return reader.ReadSingle();
case DataType.Int64:
return reader.ReadInt64();
case DataType.UInt64:
return reader.ReadUInt64();
default:
throw new NotSupportedException();
}
}
}
}
Write Preview
Loading…
Cancel
Save