text:
/* vector/gsl_vector_uchar.h
 *
 * Copyright (C) 1996, 1997, 1998, 1999, 2000 Gerard Jungman, Brian Gough
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or (at
 * your option) any later version.
 *
 * This program is distributed in the hope that it will be useful, but
 * WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
 */

#ifndef __GSL_VECTOR_UCHAR_H__
#define __GSL_VECTOR_UCHAR_H__

#include <stdlib.h>
#include <gsl/gsl_types.h>
#include <gsl/gsl_errno.h>
#include <gsl/gsl_check_range.h>
#include <gsl/gsl_block_uchar.h>

#undef __BEGIN_DECLS
#undef __END_DECLS
#ifdef __cplusplus
# define __BEGIN_DECLS extern "C" {
# define __END_DECLS }
#else
# define __BEGIN_DECLS /* empty */
# define __END_DECLS /* empty */
#endif

__BEGIN_DECLS

typedef struct
{
  size_t size;
  size_t stride;
  unsigned char *data;
  gsl_block_uchar *block;
  int owner;
}
gsl_vector_uchar;

typedef struct
{
  gsl_vector_uchar vector;
} _gsl_vector_uchar_view;

typedef _gsl_vector_uchar_view gsl_vector_uchar_view;

typedef struct
{
  gsl_vector_uchar vector;
} _gsl_vector_uchar_const_view;

typedef const _gsl_vector_uchar_const_view gsl_vector_uchar_const_view;

/* Allocation */

gsl_vector_uchar *gsl_vector_uchar_alloc (const size_t n);
gsl_vector_uchar *gsl_vector_uchar_calloc (const size_t n);

gsl_vector_uchar *gsl_vector_uchar_alloc_from_block (gsl_block_uchar * b,
                                                     const size_t offset,
                                                     const size_t n,
                                                     const size_t stride);

gsl_vector_uchar *gsl_vector_uchar_alloc_from_vector (gsl_vector_uchar * v,
                                                      const size_t offset,
                                                      const size_t n,
                                                      const size_t stride);

void gsl_vector_uchar_free (gsl_vector_uchar * v);

/* Views */

_gsl_vector_uchar_view
gsl_vector_uchar_view_array (unsigned char *v, size_t n);

_gsl_vector_uchar_view
gsl_vector_uchar_view_array_with_stride (unsigned char *base,
                                         size_t stride,
                                         size_t n);

_gsl_vector_uchar_const_view
gsl_vector_uchar_const_view_array (const unsigned char *v, size_t n);

_gsl_vector_uchar_const_view
gsl_vector_uchar_const_view_array_with_stride (const unsigned char *base,
                                               size_t stride,
                                               size_t n);

_gsl_vector_uchar_view
gsl_vector_uchar_subvector (gsl_vector_uchar *v,
                            size_t i,
                            size_t n);

_gsl_vector_uchar_view
gsl_vector_uchar_subvector_with_stride (gsl_vector_uchar *v,
                                        size_t i,
                                        size_t stride,
                                        size_t n);

_gsl_vector_uchar_const_view
gsl_vector_uchar_const_subvector (const gsl_vector_uchar *v,
                                  size_t i,
                                  size_t n);

_gsl_vector_uchar_const_view
gsl_vector_uchar_const_subvector_with_stride (const gsl_vector_uchar *v,
                                              size_t i,
                                              size_t stride,
                                              size_t n);

/* Operations */

unsigned char gsl_vector_uchar_get (const gsl_vector_uchar * v, const size_t i);
void gsl_vector_uchar_set (gsl_vector_uchar * v, const size_t i, unsigned char x);

unsigned char *gsl_vector_uchar_ptr (gsl_vector_uchar * v, const size_t i);
const unsigned char *gsl_vector_uchar_const_ptr (const gsl_vector_uchar * v, const size_t i);

void gsl_vector_uchar_set_zero (gsl_vector_uchar * v);
void gsl_vector_uchar_set_all (gsl_vector_uchar * v, unsigned char x);
int gsl_vector_uchar_set_basis (gsl_vector_uchar * v, size_t i);

int gsl_vector_uchar_fread (FILE * stream, gsl_vector_uchar * v);
int gsl_vector_uchar_fwrite (FILE * stream, const gsl_vector_uchar * v);
int gsl_vector_uchar_fscanf (FILE * stream, gsl_vector_uchar * v);
int gsl_vector_uchar_fprintf (FILE * stream, const gsl_vector_uchar * v,
                              const char *format);

int gsl_vector_uchar_memcpy (gsl_vector_uchar * dest, const gsl_vector_uchar * src);

int gsl_vector_uchar_reverse (gsl_vector_uchar * v);

int gsl_vector_uchar_swap (gsl_vector_uchar * v, gsl_vector_uchar * w);
int gsl_vector_uchar_swap_elements (gsl_vector_uchar * v, const size_t i, const size_t j);

unsigned char gsl_vector_uchar_max (const gsl_vector_uchar * v);
unsigned char gsl_vector_uchar_min (const gsl_vector_uchar * v);
void gsl_vector_uchar_minmax (const gsl_vector_uchar * v, unsigned char * min_out, unsigned char * max_out);

size_t gsl_vector_uchar_max_index (const gsl_vector_uchar * v);
size_t gsl_vector_uchar_min_index (const gsl_vector_uchar * v);
void gsl_vector_uchar_minmax_index (const gsl_vector_uchar * v, size_t * imin, size_t * imax);

int gsl_vector_uchar_add (gsl_vector_uchar * a, const gsl_vector_uchar * b);
int gsl_vector_uchar_sub (gsl_vector_uchar * a, const gsl_vector_uchar * b);
int gsl_vector_uchar_mul (gsl_vector_uchar * a, const gsl_vector_uchar * b);
int gsl_vector_uchar_div (gsl_vector_uchar * a, const gsl_vector_uchar * b);
int gsl_vector_uchar_scale (gsl_vector_uchar * a, const double x);
int gsl_vector_uchar_add_constant (gsl_vector_uchar * a, const double x);

int gsl_vector_uchar_isnull (const gsl_vector_uchar * v);

#ifdef HAVE_INLINE

extern inline
unsigned char
gsl_vector_uchar_get (const gsl_vector_uchar * v, const size_t i)
{
#if GSL_RANGE_CHECK
  if (i >= v->size)
    {
      GSL_ERROR_VAL ("index out of range", GSL_EINVAL, 0);
    }
#endif
  return v->data[i * v->stride];
}

extern inline
void
gsl_vector_uchar_set (gsl_vector_uchar * v, const size_t i, unsigned char x)
{
#if GSL_RANGE_CHECK
  if (i >= v->size)
    {
      GSL_ERROR_VOID ("index out of range", GSL_EINVAL);
    }
#endif
  v->data[i * v->stride] = x;
}

extern inline
unsigned char *
gsl_vector_uchar_ptr (gsl_vector_uchar * v, const size_t i)
{
#if GSL_RANGE_CHECK
  if (i >= v->size)
    {
      GSL_ERROR_NULL ("index out of range", GSL_EINVAL);
    }
#endif
  return (unsigned char *) (v->data + i * v->stride);
}

extern inline
const unsigned char *
gsl_vector_uchar_const_ptr (const gsl_vector_uchar * v, const size_t i)
{
#if GSL_RANGE_CHECK
  if (i >= v->size)
    {
      GSL_ERROR_NULL ("index out of range", GSL_EINVAL);
    }
#endif
  return (const unsigned char *) (v->data + i * v->stride);
}

#endif /* HAVE_INLINE */

__END_DECLS

#endif /* __GSL_VECTOR_UCHAR_H__ */

ext: h


text:
/*

   Copyright 2016  Amazon.com, Inc. or its affiliates. All Rights Reserved.

   Licensed under the Apache License, Version 2.0 (the "License"). You may not use this file except in compliance with the License. A copy of the License is located at

   http://aws.amazon.com/apache2.0/

   or in the "license" file accompanying this file. This file is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License.
 */

#ifndef NNTYPES_H
#define NNTYPES_H
#include <vector>
#include <set>
#include <string>
#include <iostream>
#include <iomanip>
#include <fstream>
#include <algorithm>
#include <netcdf>
#ifndef __NVCC__
#include <tuple>
#include <json/json.h>
#endif
#include <cmath>
#include <memory>

class NNDataSetBase;
class NNLayer;
class NNNetwork;
class NNWeight;

// Activates step by step CPU validation
#define VALIDATION
#ifdef VALIDATION
extern "C"
{
    #include <cblas.h>
}
#endif

static const float NN_VERSION       = 0.9f;
static const float MIN_ERROR        = 1.0e-12f;
static const float MIN_ACTIVATION   = 0.000001f;
static const float MAX_ACTIVATION   = 0.999999f;
static const float MAX_VALUE        = 999999999999999.0f;

template <typename T> struct GpuBuffer;

enum
{
    DefaultBatch    = 512
};

enum Mode {
    Prediction = 0,
    Training = 1,
    Validation = 2,
    Unspecified = 3
};

enum TrainingMode
{
    SGD = 0,
    Momentum = 1,
    AdaGrad = 2,
    Nesterov = 3,
    RMSProp = 4,
    AdaDelta = 5,
    Adam = 6,
};

ostream& operator<< (ostream& out, const TrainingMode& e);

enum ErrorFunction
{
    L1,
    L2,
    CrossEntropy,
    ScaledMarginalCrossEntropy,
    DataScaledMarginalCrossEntropy,
    Hinge,
    L2Hinge,
};

ostream& operator<< (ostream& out, const ErrorFunction& e);

enum Activation {
    Sigmoid,
    Tanh,
    RectifiedLinear,
    Linear,
    ParametricRectifiedLinear,
    SoftPlus,
    SoftSign,
    SoftMax,
    RELUMax,
    LinearMax,
    ExponentialLinear,
    LeakyRectifiedLinear,
    ScaledExponentialLinear,
};

ostream& operator<< (ostream& out, const Activation& a);

enum WeightInitialization
{
    Xavier,
    CaffeXavier,
    Gaussian,
    Uniform,
    UnitBall,
    Constant,
    SELU,
};

ostream& operator<< (ostream& out, const WeightInitialization& w);

enum PoolingFunction {
    None,
    Max,
    Average,
    LRN,
    Maxout,
    DotProduct,
    Cosine,
    Stochastic,
    LCN,
    GlobalTemporal,
};

ostream& operator<< (ostream& out, const PoolingFunction& p);

#include "kernels.h"
#include "GpuSort.h"
#include "NNEnum.h"
#include "NNWeight.h"
#include "NNLayer.h"
#include "NNNetwork.h"

int MPI_Bcast_string(string& s);

struct NNDataSetDimensions
{
    uint32_t _dimensions;
    uint32_t _width;
    uint32_t _height;
    uint32_t _length;

    NNDataSetDimensions();

    NNDataSetDimensions(uint32_t width, uint32_t height = 1, uint32_t length = 1);
};

struct NNDataSetDescriptor
{
    string _name;                       // dataset name
    NNDataSetEnums::DataType _dataType; // dataset type
    uint32_t _attributes;               // dataset attributes
    NNDataSetDimensions _dim;           // dataset dimensions
    uint32_t _examples;                 // number of examples in this dataset
    float _sparseDensity;               // sparseness density of this dataset

    /**
     * Only support vanilla sparse and dense datasets for now
     */
    static bool isSupported(uint32_t attributes)
    {
        using NNDataSetEnums::Attributes;

        static const vector<Attributes> SUPPORTED_ATTRIBUTES(Attributes::Sparse);
        for (auto mask : SUPPORTED_ATTRIBUTES)
        {
            if (attributes & mask)
            {
                attributes -= mask;
            }
        }
        return attributes == 0;
    }
};

/**
 * Creates a new NNDatSetBase object given the descriptor.
 * Only attributes specified by the NNDataSetDescriptor::isSupported
 * is valid. If any other attributes are specified, a std::runtime_error
 * is thrown. It is the caller's responsibility to delete the returned object
 * (not wrapping the return value with a unique_ptr to be consistent with
 * the rest of the code base).
 */
NNDataSetBase* createNNDataSet(const NNDataSetDescriptor &descriptor);

struct NNDataSetBase {

    string                          _name;                          // Dataset name
    NNDataSetEnums::DataType        _dataType;                      // Dataset type (see above enum)
    uint32_t                        _attributes;                    // Dataset characteristics (see NNDataSetEnum::Attributes in NNEnum.h)
    uint32_t                        _examples;                      // Number of examples
    uint32_t                        _uniqueExamples;                // Number of unique examples for indexed data (set to examples if unindexed)
    uint32_t                        _localExamples;                 // Number of local examples when data sharded
    uint32_t                        _dimensions;                    // Dimensionality of data set
    uint32_t                        _width;                         // Dataset x dimension
    uint32_t                        _height;                        // Dataset y dimension
    uint32_t                        _length;                        // Dataset z dimension
    uint32_t                        _stride;                        // Stride between examples
    NNDataSetEnums::Sharding        _sharding;                      // Sharding of dataset for parallel execution
    uint32_t                        _minX;                          // Beginning of local X sharding for model parallel execution
    uint32_t                        _maxX;                          // End of local X sharding for model parallel execution
    uint64_t                        _sparseDataSize;                // Total sparse datapoints
    NNFloat                         _sparseDensity;                 // Overall sparse density (0.0 - 1.0)
    vector<uint64_t>                _vSparseStart;                  // Vector of sparse datapoint starts per example
    unique_ptr<GpuBuffer<uint64_t>> _pbSparseStart;                 // GPU copy of _vSparseStart
    vector<uint64_t>                _vSparseEnd;                    // Vector of sparse datapoint ends per example
    unique_ptr<GpuBuffer<uint64_t>> _pbSparseEnd;                   // GPU copy of _vSparseEnd
    vector<uint32_t>                _vSparseIndex;                  // Vector of sparse indices
    unique_ptr<GpuBuffer<uint32_t>> _pbSparseIndex;                 // GPU copy of _vSparseIndex
    vector<NNFloat>                 _vDataWeight;                   // Per example sparse index weights
    unique_ptr<GpuBuffer<NNFloat>>  _pbDataWeight;                  // GPU copy of _vDataWeight
    vector<uint32_t>                _vIndex;                        // Indexed data array
    unique_ptr<GpuBuffer<uint32_t>> _pbIndex;                       // GPU copy of _vIndex;
    unique_ptr<GpuBuffer<NNFloat>>  _pbDenoisingRandom;             // Denoising randoms

    // Transposed sparse lookup for sparse backpropagation
    vector<uint64_t>                _vSparseDatapointCount;
    vector<uint32_t>                _vSparseMaxDatapointCount;
    vector<uint32_t>                _vSparseMultiDatapointCount;
    vector<uint32_t>                _vSparseTransposedStart;
    uint64_t                        _sparseTransposedIndices;
    unique_ptr<GpuBuffer<uint32_t>> _pbSparseTransposedStart;
    unique_ptr<GpuBuffer<uint32_t>> _pbSparseTransposedEnd;
    unique_ptr<GpuBuffer<uint32_t>> _pbSparseTransposedIndex;
    unique_ptr<GpuBuffer<NNFloat>>  _pbSparseTransposedData;

    // States
    bool                            _bDenoising;
    bool                            _bDirty;
    bool                            _bStreaming;
    bool                            _bIndexed;
    uint32_t                        _batch;

    NNDataSetBase();
    NNDataSetDimensions GetDimensions();
    uint32_t GetExamples() { return _examples; };
    uint32_t GetUniqueExamples() { return _uniqueExamples; };

    virtual bool SaveNetCDF(const string& fname) = 0;
    virtual bool WriteNetCDF(netCDF::NcFile& nfc, const string& fname, const uint32_t n) = 0;
    virtual ~NNDataSetBase() = 0;
    virtual void RefreshState(uint32_t batch) = 0;
    virtual bool Shard(NNDataSetEnums::Sharding sharding) = 0;
    virtual bool UnShard() = 0;
    virtual bool SetStreaming(bool flag) = 0;
    virtual bool GetStreaming() = 0;
    virtual vector<tuple<uint64_t, uint64_t> > getMemoryUsage() = 0;
    virtual bool CalculateSparseDatapointCounts() = 0;
    virtual bool GenerateSparseTransposedMatrix(uint32_t batch, NNLayer* pLayer) = 0;
    virtual bool CalculateSparseTransposedMatrix(uint32_t position, uint32_t batch, NNLayer* pLayer) = 0;
    virtual bool CalculateSparseTransposedDenoisedMatrix(uint32_t position, uint32_t batch, NNLayer* pLayer) = 0;
    virtual bool CalculateSparseTransposedWeightGradient(NNFloat alpha, NNFloat beta, uint32_t m, uint32_t n, NNFloat* pDelta, NNFloat* pWeightGradient) = 0;
    virtual bool SetDenoising(bool flag) = 0;
    virtual bool GenerateDenoisingData() = 0;
    virtual bool LoadInputUnit(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual bool LoadSparseInputUnit(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual bool LoadSparseDenoisedInputUnit(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual bool CalculateSparseZ(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pWeight, NNFloat* pUnit, NNFloat beta = (NNFloat)0.0) = 0;
    virtual bool CalculateSparseDenoisedZ(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pWeight, NNFloat* pUnit, NNFloat beta = (NNFloat)0.0) = 0;
    virtual float CalculateL1Error(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual float CalculateL2Error(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual float CalculateL2HingeError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual float CalculateCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual float CalculateScaledMarginalCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual float CalculateMultinomialCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual float CalculateMultinomialScaledMarginalCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual float CalculateDataScaledMarginalCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual float CalculateHingeError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit) = 0;
    virtual bool CalculateL1OutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta, NNFloat slope, NNFloat alpha, NNFloat lambda) = 0;
    virtual bool CalculateCrossEntropyOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta) = 0;
    virtual bool CalculateScaledMarginalCrossEntropyOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta) = 0;
    virtual bool CalculateOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta, NNFloat slope, NNFloat alpha, NNFloat lambda) = 0;
    virtual bool CalculateL2HingeOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta, NNFloat slope, NNFloat alpha, NNFloat lambda) = 0;
    virtual bool CalculateDataScaledMarginalCrossEntropyOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta) = 0;
    virtual bool CalculateHingeOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta) = 0;

    /**
     * Copies data from srcData to this dataset and uploads to GPU.
     * The length to copy is determined by the example size and stride
     * of this dataset. Only valid for dense and dense indexed datasets.
     * If the dataset is indexed this sets the unique data
     * and NNDataSet::SetIndexedData sets the actual examples.
     * An exception is thrown if this dataset is not dense.
     */
    virtual void LoadDenseData(const void *srcData) = 0;

    /**
     * Copies data from srcData to this dataset but doesn't upload to GPU.
     */
    virtual void CopyDenseData(const void *srcData) = 0;

    /**
     * Copies sparse data points from the specified src to this NNDataSet
     * and uploads to GPU. Any existing data is overwritten.
     * An exception is thrown if this dataset is not sparse.
     */
    virtual void LoadSparseData(const uint64_t *srcSparseStart, const uint64_t *srcSparseEnd, const void *srcSparseData,
                                const uint32_t *srcSparseIndex) = 0;

    /**
     * Copies sparse data points from the specified src to this NNDataSet
     * but doesn't upload to GPU. Any existing data is overwritten.
     * An exception is thrown if this dataset is not sparse.
     */
    virtual void CopySparseData(const uint64_t *srcSparseStart, const uint64_t *srcSparseEnd, const void *srcSparseData,
                                const uint32_t *srcSparseIndex) = 0;

    /**
     * Same as LoadSparseData(const uint64_t*, const uint64_t*, const void*, const uint32_t*) except that
     * it takes long* for the sparse start, end and indexes and casts them into
     * uint64_t and uint32_t during load. It is up to the caller to ensure that the cast is bounds safe
     * (e.g. no negative indexes). This method is useful when writing language extensions for those languages
     * that do not have unsigned primitive types (e.g. Java).
     */
    virtual void LoadSparseData(const long *srcSparseStart, const long *srcSparseEnd, const void *srcSparseData,
                                const long *srcSparseIndex) = 0;

      /**
     * Same as CopySparseData(const uint64_t*, const uint64_t*, const void*, const uint32_t*) except that
     * it takes long* for the sparse start, end and indexes and casts them into
     * uint64_t and uint32_t during load. It is up to the caller to ensure that the cast is bounds safe
     * (e.g. no negative indexes). This method is useful when writing language extensions for those languages
     * that do not have unsigned primitive types (e.g. Java).
     */
    virtual void CopySparseData(const long *srcSparseStart, const long *srcSparseEnd, const void *srcSparseData,
                                const long *srcSparseIndex) = 0;

    /**
     * If this dataset is indexed, then sets the actual (unique) examples.
     * An exception is thrown if this dataset is not indexed.
     */
    virtual void LoadIndexedData(const uint32_t *srcIndexedData) = 0;

    /**
     * If this dataset is weighted, then sets the weights for each example.
     * An exception is thrown if this dataset is not weighted.
     */
    virtual void LoadDataWeight(const NNFloat *srcWeightData) = 0;

 protected:
    NNDataSetBase(const string &name, NNDataSetEnums::DataType dataType, uint32_t examples, uint32_t uniqueExamples,
                  const NNDataSetDimensions &datasetDim);

};

ostream& operator<< (ostream& out, NNDataSetEnums::Attributes& a);
ostream& operator<< (ostream& out, NNDataSetEnums::Kind& k);
ostream& operator<< (ostream& out, NNDataSetEnums::DataType& t);
ostream& operator<< (ostream& out, NNDataSetEnums::Sharding& s);

template<typename T> class NNDataSet : public NNDataSetBase {
public:
    friend class NNetwork;
    friend class NNLayer;
    friend vector<NNDataSetBase*> LoadNetCDF(const string& fname);
    friend bool SaveNetCDF(const string& fname, vector<NNDataSetBase*> vDataSet);

private:
    // Type-specific data
    vector<T>                   _vData;
    unique_ptr<GpuBuffer<T>>    _pbData;
    vector<T>                   _vSparseData;
    unique_ptr<GpuBuffer<T>>    _pbSparseData;

    // Force constructor private
    NNDataSet(const string& fname, uint32_t n);
    bool Rename(const string& name);
    bool SaveNetCDF(const string& fname);
    bool WriteNetCDF(netCDF::NcFile& nfc, const string& fname, const uint32_t n);
    void RefreshState(uint32_t batch) {}
    bool Shard(NNDataSetEnums::Sharding sharding);
    bool UnShard();
    vector<tuple<uint64_t, uint64_t> > getMemoryUsage();
    bool CalculateSparseDatapointCounts();
    bool GenerateSparseTransposedMatrix(uint32_t batch, NNLayer* pLayer);
    bool CalculateSparseTransposedMatrix(uint32_t position, uint32_t batch, NNLayer* pLayer);
    bool CalculateSparseTransposedDenoisedMatrix(uint32_t position, uint32_t batch, NNLayer* pLayer);
    bool CalculateSparseTransposedWeightGradient(NNFloat alpha, NNFloat beta, uint32_t m, uint32_t n, NNFloat* pDelta, NNFloat* pWeightGradient);
    bool SetStreaming(bool flag);
    bool GetStreaming();
    bool SetDenoising(bool flag);
    bool GenerateDenoisingData();
    bool LoadInputUnit(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    bool LoadSparseInputUnit(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    bool LoadSparseDenoisedInputUnit(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    bool CalculateSparseZ(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pWeight, NNFloat* pUnit, NNFloat beta);
    bool CalculateSparseDenoisedZ(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pWeight, NNFloat* pUnit, NNFloat beta);
    float CalculateL1Error(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    float CalculateL2Error(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    float CalculateL2HingeError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    float CalculateCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    float CalculateScaledMarginalCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    float CalculateMultinomialCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    float CalculateMultinomialScaledMarginalCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    float CalculateDataScaledMarginalCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    float CalculateHingeError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit);
    bool CalculateL1OutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta, NNFloat slope, NNFloat alpha, NNFloat lambda);
    bool CalculateCrossEntropyOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta);
    bool CalculateScaledMarginalCrossEntropyOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta);
    bool CalculateOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta, NNFloat slope, NNFloat alpha, NNFloat lambda);
    bool CalculateL2HingeOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta, NNFloat slope, NNFloat alpha, NNFloat lambda);
    bool CalculateDataScaledMarginalCrossEntropyOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta);
    bool CalculateHingeOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta);

public:
    /**
     * Creates dense dataset with the specified dimensions and name
     * with space allocated for the specified number of examples
     */
    NNDataSet(uint32_t examples, const NNDataSetDimensions &dim, const string &name = "");

    /**
     * Creates dense indexed dataset with the specified dimensions and name
     * with space allocated for the specified number of examples
     * and index data space allocated for the specified number of uniqueExamples.
     */
    NNDataSet(uint32_t examples, uint32_t uniqueExamples, const NNDataSetDimensions &dim, const string &name = "");

    /**
     * Creates sparse dataset for the layer with space allocated
     * for the specified number of examples. The isWeighted parameter
     * can be set to true to create a sparse weighted dataset.
     */
    NNDataSet(uint32_t examples, NNFloat sparseDensity, const NNDataSetDimensions &dim, bool isWeighted = false, const string &name = "");

    /**
     * Creates a sparse dataset with the specified dimensions
     * and name with space allocated for the specified number of
     * examples and index data (if isIndexed).
     * The isIndexed and isWeighted parameters can be set to
     * create a sparse indexed, weighted dataset.
     */
    NNDataSet(uint32_t examples, uint32_t uniqueExamples, size_t sparseDataSize, const NNDataSetDimensions &dim,
              bool isIndexed = false, bool isWeighted = false, const string &name = "");

    void LoadDenseData(const void *srcData) override;

    void CopyDenseData(const void *srcData) override;

    void LoadSparseData(const uint64_t *srcSparseStart, const uint64_t *srcSparseEnd, const void *srcSparseData,
                        const uint32_t *srcSparseIndex) override;

    void CopySparseData(const uint64_t *srcSparseStart, const uint64_t *srcSparseEnd, const void *srcSparseData,
                        const uint32_t *srcSparseIndex) override;

    void LoadSparseData(const long *srcSparseStart, const long *srcSparseEnd, const void *srcSparseData,
                        const long *srcSparseIndex) override;

    void CopySparseData(const long *srcSparseStart, const long *srcSparseEnd, const void *srcSparseData,
                        const long *srcSparseIndex) override;

    void LoadIndexedData(const uint32_t *srcIndexedData) override;

    void LoadDataWeight(const NNFloat *srcWeightData) override;

    ~NNDataSet();
    void Shuffle();
    T GetDataPoint(uint32_t n, uint32_t x, uint32_t y = 0, uint32_t z = 0);
    bool SetDataPoint(T v, uint32_t n, uint32_t x, uint32_t y = 0, uint32_t z = 0);
    uint64_t GetSparseDataPoints(uint32_t n);
    uint32_t GetSparseIndex(uint32_t n, uint32_t i);
    bool SetSparseIndex(uint32_t n, uint32_t i, uint32_t v);
    T GetSparseDataPoint(uint32_t n, uint32_t i);
    bool SetSparseDataPoint(uint32_t n, uint32_t i, T v);
};

template<typename T> bool NNDataSet<T>::LoadInputUnit(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    if (_attributes & NNDataSetEnums::Indexed)
        kLoadIndexedInputUnit(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbData->_pDevData);
    else
        kLoadInputUnit(position, batch, stride, pUnit, _pbData->_pDevData);
    return true;
}

template<typename T> bool NNDataSet<T>::LoadSparseInputUnit(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Boolean)
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kLoadIndexedSparseInputUnit(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight);
        else
            kLoadSparseInputUnit(position, batch, stride, pUnit, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight);
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kLoadIndexedSparseAnalogInputUnit(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData);
        else
            kLoadSparseAnalogInputUnit(position, batch, stride, pUnit, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData);
    }
    return true;
}

template<typename T> bool NNDataSet<T>::LoadSparseDenoisedInputUnit(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Boolean)
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kLoadIndexedSparseDenoisedInputUnit(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbDenoisingRandom->_pDevData);
        else
            kLoadSparseDenoisedInputUnit(position, batch, stride, pUnit, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbDenoisingRandom->_pDevData);
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kLoadIndexedSparseAnalogDenoisedInputUnit(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, _pbDenoisingRandom->_pDevData);
        else
            kLoadSparseAnalogDenoisedInputUnit(position, batch, stride, pUnit, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, _pbDenoisingRandom->_pDevData);
    }
    return true;
}

template<typename T> bool NNDataSet<T>::CalculateSparseZ(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pWeight, NNFloat* pUnit, NNFloat beta)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Boolean)
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedSparseZ(position, batch, stride, pWeight, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, pUnit, beta);
        else
            kCalculateSparseZ(position, batch, stride, pWeight, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, pUnit, beta);
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedSparseAnalogZ(position, batch, stride, pWeight, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, pUnit, beta);
        else
            kCalculateSparseAnalogZ(position, batch, stride, pWeight, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, pUnit, beta);
    }
    return true;
}

template<typename T> bool NNDataSet<T>::CalculateSparseDenoisedZ(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pWeight, NNFloat* pUnit, NNFloat beta)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Boolean)
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedSparseDenoisedZ(position, batch, stride, pWeight, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbDenoisingRandom->_pDevData, pUnit, beta);
        else
            kCalculateSparseDenoisedZ(position, batch, stride, pWeight, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbDenoisingRandom->_pDevData, pUnit, beta);
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedSparseAnalogDenoisedZ(position, batch, stride, pWeight, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, _pbDenoisingRandom->_pDevData, pUnit, beta);
        else
            kCalculateSparseAnalogDenoisedZ(position, batch, stride, pWeight, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, _pbDenoisingRandom->_pDevData, pUnit, beta);
    }
    return true;
}

template<typename T> bool NNDataSet<T>::CalculateSparseTransposedMatrix(uint32_t position, uint32_t batch, NNLayer* pLayer)
{
    // Rebuild sparse data table if dataset changed
    if (_bDirty || (batch != _batch))
    {
        GenerateSparseTransposedMatrix(batch, pLayer);
    }

    // Initialize transposed sparse offsets
    _pbSparseTransposedEnd->Copy(_pbSparseTransposedStart->_pDevData);

    // Call appropriate matrix generation kernel
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    NNFloat* pSparseTransposedData = ((_attributes & NNDataSetEnums::Weighted) || !(_attributes & NNDataSetEnums::Boolean)) ? _pbSparseTransposedData->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Boolean)
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedSparseTransposedMatrix(position, batch, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseTransposedEnd->_pDevData, _pbSparseTransposedIndex->_pDevData, pSparseTransposedData);
        else
            kCalculateSparseTransposedMatrix(position, batch, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseTransposedEnd->_pDevData, _pbSparseTransposedIndex->_pDevData, pSparseTransposedData);
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedSparseTransposedAnalogMatrix(position, batch, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, _pbSparseTransposedEnd->_pDevData, _pbSparseTransposedIndex->_pDevData, pSparseTransposedData);
        else
            kCalculateSparseTransposedAnalogMatrix(position, batch, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, _pbSparseTransposedEnd->_pDevData, _pbSparseTransposedIndex->_pDevData, pSparseTransposedData);
    }

    return true;
}

template<typename T> bool NNDataSet<T>::CalculateSparseTransposedDenoisedMatrix(uint32_t position, uint32_t batch, NNLayer* pLayer)
{

    // Rebuild sparse data table if dataset changed
    if (_bDirty || (batch != _batch))
    {
        GenerateSparseTransposedMatrix(batch, pLayer);
    }

    // Initialize transposed sparse offsets
    _pbSparseTransposedEnd->Copy(_pbSparseTransposedStart->_pDevData);

    // Call appropriate matrix generation kernel
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    NNFloat* pSparseTransposedData = ((_attributes & NNDataSetEnums::Weighted) || !(_attributes & NNDataSetEnums::Boolean)) ? _pbSparseTransposedData->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Boolean)
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedSparseTransposedDenoisedMatrix(position, batch, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbDenoisingRandom->_pDevData, _pbSparseTransposedEnd->_pDevData, _pbSparseTransposedIndex->_pDevData, pSparseTransposedData);
        else
            kCalculateSparseTransposedDenoisedMatrix(position, batch, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbDenoisingRandom->_pDevData, _pbSparseTransposedEnd->_pDevData, _pbSparseTransposedIndex->_pDevData, pSparseTransposedData);
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedSparseTransposedAnalogDenoisedMatrix(position, batch, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, _pbDenoisingRandom->_pDevData, _pbSparseTransposedEnd->_pDevData, _pbSparseTransposedIndex->_pDevData, pSparseTransposedData);
        else
            kCalculateSparseTransposedAnalogDenoisedMatrix(position, batch, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, _pbDenoisingRandom->_pDevData, _pbSparseTransposedEnd->_pDevData, _pbSparseTransposedIndex->_pDevData, pSparseTransposedData);
    }

#if 0
    vector<uint32_t> vSparseTransposedStart(53120);
    vector<uint32_t> vSparseTransposedEnd(53120);
    _pbSparseTransposedStart->Download(&vSparseTransposedStart[0]);
    _pbSparseTransposedEnd->Download(&vSparseTransposedEnd[0]);
    for (uint32_t i = 0; i < 53120; i++)
    printf("%6u %9u %9u %9u %9u\n", i, vSparseTransposedStart[i], vSparseTransposedEnd[i], vSparseTransposedEnd[i] - vSparseTransposedStart[i], (uint32_t)_vSparseDatapointCount[i]);
    exit(-1);
#endif
    return true;
}

template<typename T> bool NNDataSet<T>::CalculateSparseTransposedWeightGradient(NNFloat alpha, NNFloat beta, uint32_t m, uint32_t n, NNFloat* pDelta, NNFloat* pWeightGradient)
{
    if ((_attributes & NNDataSetEnums::Boolean) && !(_attributes & NNDataSetEnums::Weighted))
        kCalculateSparseTransposedWeightGradient(alpha, beta, m, n, _pbSparseTransposedStart->_pDevData, _pbSparseTransposedEnd->_pDevData, _pbSparseTransposedIndex->_pDevData, pDelta, pWeightGradient);
    else
        kCalculateSparseTransposedAnalogWeightGradient(alpha, beta, m, n, _pbSparseTransposedStart->_pDevData, _pbSparseTransposedEnd->_pDevData, _pbSparseTransposedIndex->_pDevData, _pbSparseTransposedData->_pDevData, pDelta, pWeightGradient);
    return true;
}

template<typename T> float NNDataSet<T>::CalculateL1Error(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse)
    {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        if (_attributes & NNDataSetEnums::Boolean)
        {
           if (_attributes & NNDataSetEnums::Indexed)
           {
                return kCalculateIndexedSparseL1Error(position, batch, stride, pUnit,
                       _pbIndex->_pDevData,
                       _pbSparseStart->_pDevData,
                       _pbSparseEnd->_pDevData,
                       _pbSparseIndex->_pDevData,
                       pDataWeight,
                       bSparseIgnoreZero);
           }
           else
           {
                return kCalculateSparseL1Error(position, batch, stride, pUnit,
                       _pbSparseStart->_pDevData,
                       _pbSparseEnd->_pDevData,
                       _pbSparseIndex->_pDevData,
                       pDataWeight,
                       bSparseIgnoreZero);
           }
        }
        else
        {
           if (_attributes & NNDataSetEnums::Indexed)
           {
                return kCalculateIndexedSparseAnalogL1Error(position, batch, stride, pUnit,
                       _pbIndex->_pDevData,
                       _pbSparseStart->_pDevData,
                       _pbSparseEnd->_pDevData,
                       _pbSparseIndex->_pDevData,
                        pDataWeight,
                       _pbSparseData->_pDevData,
                       bSparseIgnoreZero);
           }
           else
           {
                return kCalculateSparseAnalogL1Error(position, batch, stride, pUnit,
                       _pbSparseStart->_pDevData,
                       _pbSparseEnd->_pDevData,
                       _pbSparseIndex->_pDevData,
                        pDataWeight,
                       _pbSparseData->_pDevData,
                       bSparseIgnoreZero);
           }
        }
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            return kCalculateIndexedL1Error(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight);
        else
            return kCalculateL1Error(position, batch, stride, pUnit, _pbData->_pDevData, pDataWeight);
    }
}

template<typename T> float NNDataSet<T>::CalculateL2Error(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse)
    {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        if (_attributes & NNDataSetEnums::Boolean)
        {
            if (_attributes & NNDataSetEnums::Indexed)
            {
                return kCalculateIndexedSparseL2Error(position, batch, stride, pUnit,
                       _pbIndex->_pDevData,
                       _pbSparseStart->_pDevData,
                       _pbSparseEnd->_pDevData,
                       _pbSparseIndex->_pDevData,
                       pDataWeight,
                       bSparseIgnoreZero);
            }
            else
            {
                return kCalculateSparseL2Error(position, batch, stride, pUnit,
                       _pbSparseStart->_pDevData,
                       _pbSparseEnd->_pDevData,
                       _pbSparseIndex->_pDevData,
                        pDataWeight,
                       bSparseIgnoreZero);
            }
        }
        else
        {
            if (_attributes & NNDataSetEnums::Indexed)
            {
                return kCalculateIndexedSparseAnalogL2Error(position, batch, stride, pUnit,
                           _pbIndex->_pDevData,
                           _pbSparseStart->_pDevData,
                           _pbSparseEnd->_pDevData,
                           _pbSparseIndex->_pDevData,
                           pDataWeight,
                           _pbSparseData->_pDevData,
                           bSparseIgnoreZero);
            }
            else
            {
                return kCalculateSparseAnalogL2Error(position, batch, stride, pUnit,
                           _pbSparseStart->_pDevData,
                           _pbSparseEnd->_pDevData,
                           _pbSparseIndex->_pDevData,
                           pDataWeight,
                           _pbSparseData->_pDevData,
                           bSparseIgnoreZero);
            }
        }
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            return kCalculateIndexedL2Error(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight);
        else
            return kCalculateL2Error(position, batch, stride, pUnit, _pbData->_pDevData, pDataWeight);
    }
}

template<typename T> float NNDataSet<T>::CalculateL2HingeError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse)
    {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        if (_attributes & NNDataSetEnums::Boolean)
        {
            if (_attributes & NNDataSetEnums::Indexed)
            {
                return kCalculateIndexedSparseL2HingeError(position, batch, stride, pUnit,
                       _pbIndex->_pDevData,
                       _pbSparseStart->_pDevData,
                       _pbSparseEnd->_pDevData,
                       _pbSparseIndex->_pDevData,
                       pDataWeight,
                       bSparseIgnoreZero);
            }
            else
            {
                return kCalculateSparseL2HingeError(position, batch, stride, pUnit,
                       _pbSparseStart->_pDevData,
                       _pbSparseEnd->_pDevData,
                       _pbSparseIndex->_pDevData,
                        pDataWeight,
                       bSparseIgnoreZero);
            }
        }
        else
        {
            if (_attributes & NNDataSetEnums::Indexed)
            {
                return kCalculateIndexedSparseAnalogL2HingeError(position, batch, stride, pUnit,
                           _pbIndex->_pDevData,
                           _pbSparseStart->_pDevData,
                           _pbSparseEnd->_pDevData,
                           _pbSparseIndex->_pDevData,
                           pDataWeight,
                           _pbSparseData->_pDevData,
                           bSparseIgnoreZero);
            }
            else
            {
                return kCalculateSparseAnalogL2HingeError(position, batch, stride, pUnit,
                           _pbSparseStart->_pDevData,
                           _pbSparseEnd->_pDevData,
                           _pbSparseIndex->_pDevData,
                           pDataWeight,
                           _pbSparseData->_pDevData,
                           bSparseIgnoreZero);
            }
        }
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            return kCalculateIndexedL2HingeError(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight);
        else
            return kCalculateL2HingeError(position, batch, stride, pUnit, _pbData->_pDevData, pDataWeight);
    }
}

template<typename T> float NNDataSet<T>::CalculateCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse)
    {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        if (_attributes & NNDataSetEnums::Indexed)
        {
            return kCalculateIndexedSparseCrossEntropyError(position, batch, stride, pUnit,
                   _pbIndex->_pDevData,
                   _pbSparseStart->_pDevData,
                   _pbSparseEnd->_pDevData,
                   _pbSparseIndex->_pDevData,
                   pDataWeight,
                   bSparseIgnoreZero);
        }
        else
        {
            return kCalculateSparseCrossEntropyError(position, batch, stride, pUnit,
                   _pbSparseStart->_pDevData,
                   _pbSparseEnd->_pDevData,
                   _pbSparseIndex->_pDevData,
                   pDataWeight,
                   bSparseIgnoreZero);
        }
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            return kCalculateIndexedCrossEntropyError(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight);
        else
            return kCalculateCrossEntropyError(position, batch, stride, pUnit, _pbData->_pDevData, pDataWeight);
    }
}

template<typename T> float NNDataSet<T>::CalculateScaledMarginalCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse)
    {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        if (_attributes & NNDataSetEnums::Indexed)
        {
            return kCalculateIndexedSparseScaledMarginalCrossEntropyError(position, batch, stride, pUnit,
                   _pbIndex->_pDevData,
                   _pbSparseStart->_pDevData,
                   _pbSparseEnd->_pDevData,
                   _pbSparseIndex->_pDevData,
                   pDataWeight,
                   bSparseIgnoreZero);
        }
    else
        {
            return kCalculateSparseScaledMarginalCrossEntropyError(position, batch, stride, pUnit,
                    _pbSparseStart->_pDevData,
                    _pbSparseEnd->_pDevData,
                    _pbSparseIndex->_pDevData,
                    pDataWeight,
                    bSparseIgnoreZero);
        }
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            return kCalculateIndexedScaledMarginalCrossEntropyError(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight);
        else
            return kCalculateScaledMarginalCrossEntropyError(position, batch, stride, pUnit, _pbData->_pDevData, pDataWeight);
    }
}

template<typename T> float NNDataSet<T>::CalculateMultinomialCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse)
    {
        if (_attributes & NNDataSetEnums::Boolean)
        {
            if (_attributes & NNDataSetEnums::Indexed)
            {
                return kCalculateIndexedSparseMultinomialCrossEntropyError(position, batch, stride, pUnit,
                        _pbIndex->_pDevData,
                        _pbSparseStart->_pDevData,
                        _pbSparseEnd->_pDevData,
                        _pbSparseIndex->_pDevData,
                        pDataWeight);
            }
            else
            {
                return kCalculateSparseMultinomialCrossEntropyError(position, batch, stride, pUnit,
                        _pbSparseStart->_pDevData,
                        _pbSparseEnd->_pDevData,
                        _pbSparseIndex->_pDevData,
                        pDataWeight);
            }
        }
        else
        {
            if (_attributes & NNDataSetEnums::Indexed)
            {
                return kCalculateIndexedSparseAnalogMultinomialCrossEntropyError(position, batch, stride, pUnit,
                        _pbIndex->_pDevData,
                        _pbSparseStart->_pDevData,
                        _pbSparseEnd->_pDevData,
                        _pbSparseIndex->_pDevData,
                        pDataWeight,
                        _pbSparseData->_pDevData);
            }
            else
            {
                return kCalculateSparseAnalogMultinomialCrossEntropyError(position, batch, stride, pUnit,
                        _pbSparseStart->_pDevData,
                        _pbSparseEnd->_pDevData,
                        _pbSparseIndex->_pDevData,
                        pDataWeight,
                        _pbSparseData->_pDevData);
            }
        }
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            return kCalculateIndexedMultinomialCrossEntropyError(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight);
        else
            return kCalculateMultinomialCrossEntropyError(position, batch, stride, pUnit, _pbData->_pDevData, pDataWeight);
    }
}

template<typename T> float NNDataSet<T>::CalculateMultinomialScaledMarginalCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse)
    {
        if (_attributes & NNDataSetEnums::Boolean)
        {
            if (_attributes & NNDataSetEnums::Indexed)
            {
                return kCalculateIndexedSparseMultinomialScaledMarginalCrossEntropyError(position, batch, stride, pUnit,
                        _pbIndex->_pDevData,
                        _pbSparseStart->_pDevData,
                        _pbSparseEnd->_pDevData,
                        _pbSparseIndex->_pDevData,
                        pDataWeight);
            }
            else
            {
                return kCalculateSparseMultinomialScaledMarginalCrossEntropyError(position, batch, stride, pUnit,
                        _pbSparseStart->_pDevData,
                        _pbSparseEnd->_pDevData,
                        _pbSparseIndex->_pDevData,
                        pDataWeight);
            }
        }
        else
        {
            if (_attributes & NNDataSetEnums::Indexed)
            {
                return kCalculateIndexedSparseAnalogMultinomialScaledMarginalCrossEntropyError(position, batch, stride, pUnit,
                        _pbIndex->_pDevData,
                        _pbSparseStart->_pDevData,
                        _pbSparseEnd->_pDevData,
                        _pbSparseIndex->_pDevData,
                        pDataWeight,
                        _pbSparseData->_pDevData);
            }
            else
            {
                return kCalculateSparseAnalogMultinomialScaledMarginalCrossEntropyError(position, batch, stride, pUnit,
                        _pbSparseStart->_pDevData,
                        _pbSparseEnd->_pDevData,
                        _pbSparseIndex->_pDevData,
                        pDataWeight,
                        _pbSparseData->_pDevData);
            }
        }
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            return kCalculateIndexedMultinomialScaledMarginalCrossEntropyError(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight);
        else
            return kCalculateMultinomialScaledMarginalCrossEntropyError(position, batch, stride, pUnit, _pbData->_pDevData, pDataWeight);
    }
}

template<typename T> float NNDataSet<T>::CalculateDataScaledMarginalCrossEntropyError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse)
    {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        // Scale by 1 if data is Boolean
        if (_attributes & NNDataSetEnums::Boolean)
        {
            if (_attributes & NNDataSetEnums::Indexed)
            {
                return kCalculateIndexedSparseScaledMarginalCrossEntropyError(position, batch, stride, pUnit,
                        _pbIndex->_pDevData,
                        _pbSparseStart->_pDevData,
                        _pbSparseEnd->_pDevData,
                        _pbSparseIndex->_pDevData,
                        pDataWeight,
                        bSparseIgnoreZero);
            }
            else
            {
                return kCalculateSparseScaledMarginalCrossEntropyError(position, batch, stride, pUnit,
                        _pbSparseStart->_pDevData,
                        _pbSparseEnd->_pDevData,
                        _pbSparseIndex->_pDevData,
                        pDataWeight,
                        bSparseIgnoreZero);
            }
        }
        else
        {
            if (_attributes & NNDataSetEnums::Indexed)
            {
                return kCalculateIndexedSparseDataScaledMarginalCrossEntropyError(position, batch, stride, pUnit,
                                _pbIndex->_pDevData,
                                _pbSparseStart->_pDevData,
                                _pbSparseEnd->_pDevData,
                                _pbSparseIndex->_pDevData,
                                _pbSparseData->_pDevData, bSparseIgnoreZero);
            }
            else
            {
                return kCalculateSparseDataScaledMarginalCrossEntropyError(position, batch, stride, pUnit,
                                _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData,
                                _pbSparseData->_pDevData, bSparseIgnoreZero);
            }
        }
    }
    else
    {
        // BUG BUG BUG unacceptable behavior, this should be caught at startup
        cout << "unsupported data format of this cost function" << endl;
        getGpu().Shutdown();
        exit(-1);
    }
}

template<typename T> float NNDataSet<T>::CalculateHingeError(uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Indexed)
        return kCalculateIndexedHingeError(position, batch, stride, pUnit, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight);
    else
        return kCalculateHingeError(position, batch, stride, pUnit, _pbData->_pDevData, pDataWeight);
}

template<typename T> bool NNDataSet<T>::CalculateL1OutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta, NNFloat slope, NNFloat alpha, NNFloat lambda)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse)
    {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedSparseL1OutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, bSparseIgnoreZero, slope, alpha, lambda);
        else
            kCalculateSparseL1OutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, bSparseIgnoreZero, slope, alpha, lambda);
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedL1OutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight, slope, alpha, lambda);
        else
            kCalculateL1OutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbData->_pDevData, pDataWeight, slope, alpha, lambda);
    }
    return true;
}

template<typename T> bool NNDataSet<T>::CalculateCrossEntropyOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse)
    {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedSparseCrossEntropyOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, bSparseIgnoreZero);
        else
            kCalculateSparseCrossEntropyOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, bSparseIgnoreZero);

    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedCrossEntropyOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight);
        else
            kCalculateCrossEntropyOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbData->_pDevData, pDataWeight);
    }
    return true;
}

template<typename T> bool NNDataSet<T>::CalculateScaledMarginalCrossEntropyOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse)
    {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedSparseScaledMarginalCrossEntropyOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, bSparseIgnoreZero);
        else
            kCalculateSparseScaledMarginalCrossEntropyOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, bSparseIgnoreZero);
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedScaledMarginalCrossEntropyOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight);
        else
            kCalculateScaledMarginalCrossEntropyOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbData->_pDevData, pDataWeight);
    }
    return true;
}

template<typename T> bool NNDataSet<T>::CalculateOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta, NNFloat slope, NNFloat alpha, NNFloat lambda)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse) {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        if (_attributes & NNDataSetEnums::Boolean)
        {
            if (_attributes & NNDataSetEnums::Indexed)
                kCalculateIndexedSparseOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, bSparseIgnoreZero, slope, alpha, lambda);
            else
                kCalculateSparseOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, bSparseIgnoreZero, slope, alpha, lambda);
        }
        else
        {
            if (_attributes & NNDataSetEnums::Indexed)
                kCalculateIndexedSparseAnalogOutputDelta(activation, position, batch, stride, pUnit,  pDelta, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, bSparseIgnoreZero, slope, alpha, lambda);
            else
                kCalculateSparseAnalogOutputDelta(activation, position, batch, stride, pUnit,  pDelta, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, bSparseIgnoreZero, slope, alpha, lambda);
        }
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight, slope, alpha, lambda);
        else
            kCalculateOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbData->_pDevData, pDataWeight, slope, alpha, lambda);
    }
    return true;
}

template<typename T> bool NNDataSet<T>::CalculateL2HingeOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta, NNFloat slope, NNFloat alpha, NNFloat lambda)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Sparse) {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        if (_attributes & NNDataSetEnums::Boolean)
        {
            if (_attributes & NNDataSetEnums::Indexed)
                kCalculateIndexedSparseL2HingeOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, bSparseIgnoreZero, slope, alpha, lambda);
            else
                kCalculateSparseL2HingeOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, bSparseIgnoreZero, slope, alpha, lambda);
        }
        else
        {
            if (_attributes & NNDataSetEnums::Indexed)
                kCalculateIndexedSparseAnalogL2HingeOutputDelta(activation, position, batch, stride, pUnit,  pDelta, _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, bSparseIgnoreZero, slope, alpha, lambda);
            else
                kCalculateSparseAnalogL2HingeOutputDelta(activation, position, batch, stride, pUnit,  pDelta, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData, pDataWeight, _pbSparseData->_pDevData, bSparseIgnoreZero, slope, alpha, lambda);
        }
    }
    else
    {
        if (_attributes & NNDataSetEnums::Indexed)
            kCalculateIndexedL2HingeOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight, slope, alpha, lambda);
        else
            kCalculateL2HingeOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbData->_pDevData, pDataWeight, slope, alpha, lambda);
    }
    return true;
}

template<typename T> bool NNDataSet<T>::CalculateDataScaledMarginalCrossEntropyOutputDelta(Activation activation,
                uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta)
{
    if (_attributes & NNDataSetEnums::Sparse)
    {
        bool bSparseIgnoreZero = _attributes & NNDataSetEnums::SparseIgnoreZero;
        if (_attributes & NNDataSetEnums::Indexed)
        {
            kCalculateIndexedSparseDataScaledMarginalCrossEntropyOutputDelta(activation, position, batch, stride, pUnit, pDelta,
                             _pbIndex->_pDevData, _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData,
                            _pbSparseData->_pDevData, bSparseIgnoreZero);
        }
        else
        {
            kCalculateSparseDataScaledMarginalCrossEntropyOutputDelta(activation, position, batch, stride, pUnit, pDelta,
                            _pbSparseStart->_pDevData, _pbSparseEnd->_pDevData, _pbSparseIndex->_pDevData,
                            _pbSparseData->_pDevData, bSparseIgnoreZero);
        }
    }
    else
    {
        cout << "unsupported data format of this cost function" << endl;
        getGpu().Shutdown();
        exit(-1);
    }
    return true;
}

template<typename T> bool NNDataSet<T>::CalculateHingeOutputDelta(Activation activation, uint32_t position, uint32_t batch, uint32_t stride, NNFloat* pUnit, NNFloat* pDelta)
{
    NNFloat* pDataWeight = (_attributes & NNDataSetEnums::Weighted) ? _pbDataWeight->_pDevData : NULL;
    if (_attributes & NNDataSetEnums::Indexed)
        kCalculateIndexedHingeOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbIndex->_pDevData, _pbData->_pDevData, pDataWeight);
    else
        kCalculateHingeOutputDelta(activation, position, batch, stride, pUnit, pDelta, _pbData->_pDevData, pDataWeight);
    return true;
}

vector<NNDataSetBase*> LoadNetCDF(const string& fname);
bool SaveNetCDF(const string& fname, vector<NNDataSetBase*> vDataset);
vector<NNDataSetBase*> LoadImageData(const string& fname);
vector<NNDataSetBase*> LoadCSVData(const string& fname);
vector<NNDataSetBase*> LoadJSONData(const string& fname);
vector<NNDataSetBase*> LoadAudioData(const string& name);

#endif

max_stars_repo_name: farruhnet/AmazonDeepEngine

max_issues_repo_name: saminfante/amazon-dsstne

max_forks_repo_name: saminfante/amazon-dsstne


text:
/* ----------------------------------------------------------------------------
  JAM_AXI_VEL_WMMT

    Calculates weighted first moments.

    INPUTS
      xp    : projected x' [pc]
      yp    : projected y' [pc]
      nxy   : number of x' and y' values given
      incl  : inclination [radians]
      lum   : projected luminous MGE
      pot   : projected potential MGE
      beta  : velocity anisotropy (1 - vz^2 / vr^2)
      kappa : rotation parameter

    NOTES
      * Based on janis1_weighted_first_moment IDL code by Michele Cappellari.

  Laura L Watkins [lauralwatkins@gmail.com]

  This code is released under a BSD 2-clause license.
  If you use this code for your research, please cite:
  Watkins et al. 2013, MNRAS, 436, 2598
  "Discrete dynamical models of omega Centauri"
  http://adsabs.harvard.edu/abs/2013MNRAS.436.2598W
---------------------------------------------------------------------------- */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <gsl/gsl_integration.h>
#include <gsl/gsl_errno.h>
#include "jam.h"
#include "../mge/mge.h"
#include "../tools/tools.h"

double** jam_axi_vel_wmmt( double *xp, double *yp, int nxy, double incl, \
        struct multigaussexp *lum, struct multigaussexp *pot, \
        double *beta, double *kappa, int *gslFlag_vel ) {

    struct params_losint lp;
    struct multigaussexp ilum, ipot;
    double *bani, *s2l, *q2l, *s2q2l, *s2p, *e2p;
    double *iz0, *iz1;
    double lim, result, error, si, ci, trpig, **sb_mu1;
    int i;
    size_t neval;
    int status = 0, *gslFlag;
    int TadejaVar = 0; gslFlag = &TadejaVar;

    // ---------------------------------

    // MGE parameters for inner integrand function

    // convert from projected MGEs to intrinsic MGEs
    ilum = mge_deproject( lum, incl );
    ipot = mge_deproject( pot, incl );

    // calculate useful quantities
    bani = (double *) malloc( ilum.ntotal * sizeof( double ) );
    s2l = (double *) malloc( ilum.ntotal * sizeof( double ) );
    q2l = (double *) malloc( ilum.ntotal * sizeof( double ) );
    s2q2l = (double *) malloc( ilum.ntotal * sizeof( double ) );
    s2p = (double *) malloc( ipot.ntotal * sizeof( double ) );
    e2p = (double *) malloc( ipot.ntotal * sizeof( double ) );

    for ( i = 0; i < ilum.ntotal; i++ ) {
        bani[i] = 1. / ( 1. - beta[i] );
        s2l[i] = pow( ilum.sigma[i], 2 );
        q2l[i] = pow( ilum.q[i], 2 );
        s2q2l[i] = s2l[i] * q2l[i];
    }

    for ( i = 0; i < ipot.ntotal; i++ ) {
        s2p[i] = pow( ipot.sigma[i], 2 );
        e2p[i] = 1. - pow( ipot.q[i], 2 );
    }

    // parameters for integrand function
    lp.incl = incl;
    lp.lum = &ilum;
    lp.pot = &ipot;
    lp.bani = bani;
    lp.s2l = s2l;
    lp.q2l = q2l;
    lp.s2q2l = s2q2l;
    lp.s2p = s2p;
    lp.e2p = e2p;
    lp.kappa = kappa;
    lp.gslFlag_losint = gslFlag;

    // ---------------------------------

    // set up integration
    gsl_integration_workspace *w = gsl_integration_workspace_alloc( 1000 );
    gsl_integration_cquad_workspace *ww = gsl_integration_cquad_workspace_alloc(1000);
    gsl_function F;
    F.function = &jam_axi_vel_losint;

    // trig angles
    si = sin( incl );
    ci = cos( incl );

    // outer limit of integration
    lim = 4. * maximum( ilum.sigma, ilum.ntotal );

    iz0 = (double *) malloc( nxy * sizeof( double ) );
    iz1 = (double *) malloc( nxy * sizeof( double ) );
    for ( i = 0; i < nxy; i++ ) {

        // parameters for integrand function
        lp.xp = xp[i];
        lp.yp = yp[i];

        // do z^0 integral
        lp.zpow = 0.;
        F.function = &jam_axi_vel_losint;
        F.params = &lp;
        gsl_integration_cquad(&F, -lim, lim, 0., 1e-4, ww, &result,
            &error, &neval);
        iz0[i] = result;

        // do z^1 integral
        lp.zpow = 1.;
        F.function = &jam_axi_vel_losint;
        F.params = &lp;
        gsl_set_error_handler_off();
        status += gsl_integration_qag( &F, -lim, lim, 1., 1., 1000, 2, w,
            &result, &error );
         /* uncomment these two lines for testing
        int test = 1;// Tadeja
        status += test; // Tadeja
        */
        if (status > 0 && *gslFlag<20) { // the limiting value is very arbitrary
            *gslFlag_vel += 1;   }

        if (status > 0 && *gslFlag>=20) { // the limiting value is very arbitrary so that the values do not exceed allowed size for integers
            *gslFlag_vel -= 1;   }
        *gslFlag_vel += *gslFlag; // if there was an integration error in losint

        if ( fabs( result ) > 1e-6 ) {
            gsl_integration_cquad(&F, -lim, lim, 0., 1e-3, ww, &result,
                &error, &neval);
        }
        iz1[i] = result;
    }

    // tidy up
    gsl_integration_workspace_free( w );
    gsl_integration_cquad_workspace_free(ww);

    // ---------------------------------

    sb_mu1 = (double **) malloc( nxy * sizeof( double* ) );
    for ( i = 0; i < nxy; i++ ) \
        sb_mu1[i] = (double *) malloc( 3 * sizeof( double ) );

    // calculate for each velocity component
    trpig = 2. * sqrt( M_PI * G );
    for ( i = 0; i < nxy; i++ ) {
        sb_mu1[i][0] = trpig * ( yp[i] * ci * iz0[i] - si * iz1[i] );
        sb_mu1[i][1] = -trpig * xp[i] * ci * iz0[i];
        sb_mu1[i][2] = trpig * xp[i] * si * iz0[i];
    }

    // ---------------------------------

    free( bani );
    free( s2l );
    free( q2l );
    free( s2q2l );
    free( s2p );
    free( e2p );

    free( iz0 );
    free( iz1 );

    return sb_mu1;

}

ext: c

max_stars_repo_name: tadejaversic/cjam

max_issues_repo_name: tadejaversic/cjam

max_forks_repo_name: tadejaversic/cjam


text:
#pragma once
/*
 * (C) Copyright 2020-2021 UCAR
 *
 * This software is licensed under the terms of the Apache Licence Version 2.0
 * which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
 */
/*! \defgroup ioda_cxx_attribute Attributes and Has_Attributes
 * \brief Ancillary data attached to variables and groups.
 * \ingroup ioda_cxx_api
 *
 * @{
 * \file Attribute.h
 * \brief @link ioda_cxx_attribute Interfaces @endlink for ioda::Attribute and related classes.
 */
#include <functional>
#include <gsl/gsl-lite.hpp>
#include <iostream>
#include <memory>
#include <string>
#include <valarray>
#include <vector>

#include "ioda/Exception.h"
#include "ioda/Misc/Dimensions.h"
#include "ioda/Misc/Eigen_Compat.h"
#include "ioda/Python/Att_ext.h"
#include "ioda/Types/Marshalling.h"
#include "ioda/Types/Type.h"
#include "ioda/Types/Type_Provider.h"
#include "ioda/defs.h"

namespace ioda {
class Attribute;

namespace detail {
class Attribute_Backend;
class Variable_Backend;

/// \brief Base class for Attributes
/// \ingroup ioda_cxx_attribute
///
/// \details You might wonder why we have this class
/// as a template. This is because we are using
/// a bit of compile-time template polymorphism to return
/// Attribute objects from base classes before Attribute is fully declared.
/// This is a variation of the (Curiously Recurring Template
/// Pattern)[https://en.wikipedia.org/wiki/Curiously_recurring_template_pattern].
template <class Attribute_Implementation = Attribute>
class Attribute_Base {
protected:
  /// Using an opaque object to implement the backend.
  std::shared_ptr<Attribute_Backend> backend_;

  // Variable_Backend's derived classes do occasionally need access to the Attribute_Backend object.
  friend class Variable_Backend;

  /// @name General Functions
  /// @{

  Attribute_Base(std::shared_ptr<Attribute_Backend>);

public:
  virtual ~Attribute_Base();

  /// @}

  /// @name Writing Data
  /// \note Writing metadata is an all-or-nothing-process, unlike writing
  /// segments of data to a variable.
  /// \note Dimensions are fixed. Attribute are not resizable.
  /// @{
  ///

  /// \brief The fundamental write function. Backends overload this function to implement all write
  /// operations.
  ///
  /// \details This function writes a span of bytes (characters) to the backend attribute storage.
  /// No type conversions take place here (see the templated conversion function, below).
  ///
  /// \param data is a span of data.
  /// \param in_memory_datatype is an opaque (backend-level) object that describes the placement of
  ///   the data in memory. Usually ignorable - needed for complex data structures.
  /// \throws ioda::Exception if data has the wrong size.
  /// \returns The attribute (for chaining).
  virtual Attribute_Implementation write(gsl::span<const char> data, const Type& type);

  /// \brief Write data.
  /// \tparam DataType is the type of the data. I.e. float, int, int32_t, uint16_t, std::string, etc.
  /// \tparam Marshaller is the class that serializes the data type into something that the backend
  /// library can use.
  /// \tparam TypeWrapper is a helper class that creates Type objects for the backend.
  /// \param data is a gsl::span (a pointer-length pair) that contains the data to be written.
  /// \param in_memory_dataType is the memory layout needed to parse data's type.
  /// \returns Another instance of this Attribute. Used for operation chaining.
  /// \throws ioda::Exception if data.size() does not match getDimensions().numElements.
  /// \see gsl::span for details of how to make a span.
  /// \see gsl::make_span
  template <class DataType, class Marshaller = ioda::Object_Accessor<DataType>,
            class TypeWrapper = Types::GetType_Wrapper<DataType>>
  Attribute_Implementation write(gsl::span<const DataType> data) {
    try {
      Marshaller m;
      auto d   = m.serialize(data);
      auto spn = gsl::make_span<const char>(reinterpret_cast<const char*>(d->DataPointers.data()),
                                            d->DataPointers.size() * Marshaller::bytesPerElement_);
      write(spn, TypeWrapper::GetType(getTypeProvider()));
      return Attribute_Implementation{backend_};
    } catch (...) {
      std::throw_with_nested(Exception(ioda_Here()));
    }
  }

  /// \brief Write data
  /// \note Normally the gsl::span write is fine. This one exists for easy Python binding.
  /// \tparam DataType is the type of the data. I.e. float, int, int32_t, uint16_t, std::string, etc.
  /// \tparam Marshaller is the class that serializes the data type into something that the backend
  /// library can use.
  /// \tparam TypeWrapper is a helper class that creates Type objects for the backend.
  /// \param data is a std::vector that contains the data to be written.
  /// \param in_memory_dataType is the memory layout needed to parse data's type.
  /// \returns Another instance of this Attribute. Used for operation chaining.
  /// \throws ioda::Exception if data.size() does not match getDimensions().numElements.
  /// \see gsl::span for details of how to make a span.
  /// \see gsl::make_span for details on how to make a span.
  template <class DataType, class Marshaller = ioda::Object_Accessor<DataType>,
            class TypeWrapper = Types::GetType_Wrapper<DataType>>
  Attribute_Implementation write(const std::vector<DataType>& data) {
    std::vector<DataType> vd = data;
    return this->write<DataType, Marshaller, TypeWrapper>(gsl::make_span(vd));
  }

  /// \brief Write data.
  /// \tparam DataType is the type of the data. I.e. float, int, int32_t, uint16_t, std::string, etc.
  /// \param data is an initializer list that contains the data to be written.
  /// \param in_memory_dataType is the memory layout needed to parse data's type.
  /// \returns Another instance of this Attribute. Used for operation chaining.
  /// \throws ioda::Exception if data.size() does not match getDimensions().numElements.
  /// \see gsl::span for details of how to make a span.
  template <class DataType>
  Attribute_Implementation write(std::initializer_list<DataType> data) {
    std::vector<DataType> v(data);
    return this->write<DataType>(gsl::make_span(v));
  }

  /// \brief Write a datum.
  /// \tparam DataType is the type of the data. I.e. float, int, int32_t, uint16_t, std::string, etc.
  /// \param data is the data to be written.
  /// \returns Another instance of this Attribute. Used for operation chaining.
  /// \throws ioda::Exception if the Attribute dimensions are larger than a single point.
  template <class DataType>  //, class Marshaller = HH::Types::Object_Accessor<DataType> >
  Attribute_Implementation write(DataType data) {
    try {
      if (getDimensions().numElements != 1)
        throw Exception("Wrong number of elements. Use a different write() method.", ioda_Here());
      return write<DataType>(gsl::make_span<DataType>(&data, 1));
    } catch (...) {
      std::throw_with_nested(Exception(ioda_Here()));
    }
  }

  /// \brief Write an Eigen object (a Matrix, an Array, a Block, a Map).
  /// \tparam EigenClass is the Eigen object to write.
  /// \param d is the data to be written.
  /// \throws ioda::Exception on a dimension mismatch.
  /// \returns the attribute
  template <class EigenClass>
  Attribute_Implementation writeWithEigenRegular(const EigenClass& d) {
#if 1  //__has_include("Eigen/Dense")
    try {
      typedef typename EigenClass::Scalar ScalarType;
      // If d is already in Row Major form, then this is optimized out.
      Eigen::Array<ScalarType, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> dout;
      dout.resize(d.rows(), d.cols());
      dout               = d;
      const auto& dconst = dout;  // To make some compilers happy.
      auto sp            = gsl::make_span(dconst.data(), static_cast<int>(d.rows() * d.cols()));

      return write<ScalarType>(sp);
    } catch (...) {
      std::throw_with_nested(Exception(ioda_Here()));
    }
#else
    static_assert(false, "The Eigen headers cannot be found, so this function cannot be used.");
#endif
  }

  /// \brief Write an Eigen Tensor-like object
  /// \tparam EigenClass is the Eigen tensor to write.
  /// \param d is the data to be written.
  /// \throws ioda::Exception on a dimension mismatch.
  /// \returns the attribute
  template <class EigenClass>
  Attribute_Implementation writeWithEigenTensor(const EigenClass& d) {
#if 1  //__has_include("unsupported/Eigen/CXX11/Tensor")
    try {
      ioda::Dimensions dims = detail::EigenCompat::getTensorDimensions(d);

      auto sp  = (gsl::make_span(d.data(), dims.numElements));
      auto res = write(sp);
      return res;
    } catch (...) {
      std::throw_with_nested(Exception(ioda_Here()));
    }
#else
    static_assert(
      false, "The Eigen unsupported/ headers cannot be found, so this function cannot be used.");
#endif
  }

  /// @}
  /// @name Reading Data
  /// @{
  ///

  /// \brief The fundamental read function. Backends overload this function to implement all read
  /// operations.
  ///
  /// \details This function reads in a span of characters from the backend attribute storage.
  /// No type conversions take place here (see the templated conversion function, below).
  ///
  /// \param data is a span of data that has length of getStorageSize().
  /// \param in_memory_datatype is an opaque (backend-level) object that describes the placement of
  ///   the data in memory. Usually ignorable - needed for complex data structures.
  /// \throws ioda::Exception if data has the wrong size.
  /// \returns The attribute (for chaining).
  virtual Attribute_Implementation read(gsl::span<char> data, const Type& in_memory_dataType) const;

  /// \brief Read data.
  ///
  /// \details This is a fundamental function that reads a span of characters from backend storage,
  /// and then performs the appropriate type conversion / deserialization into objects in data.
  ///
  /// \tparam DataType is the type if the data to be read. I.e. float, int, int32_t, uint16_t,
  /// std::string, etc.
  /// \tparam Marshaller is the class that performs the deserialization operation.
  /// \tparam TypeWrapper is a helper class that passes Type information to the backend.
  /// \param data is a pointer-size pair to the data buffer that is filled with the metadata's
  /// contents. It should be
  ///   pre-sized to accomodate all of the matadata. See getDimensions().numElements. data will be
  ///   filled in row-major order.
  /// \param in_memory_datatype is an opaque (backend-level) object that describes the placement of
  ///   the data in memory. Usually this does not need to be set to anything other than its default
  ///   value. Kept as a parameter for debugging purposes.
  /// \returns Another instance of this Attribute. Used for operation chaining.
  /// \throws ioda::Exception if data.size() != getDimensions().numElements.
  /// \see getDimensions for buffer size information.
  template <class DataType, class Marshaller = ioda::Object_Accessor<DataType>,
            class TypeWrapper = Types::GetType_Wrapper<DataType>>
  Attribute_Implementation read(gsl::span<DataType> data) const {
    try {
      const size_t numObjects = data.size();
      if (getDimensions().numElements != gsl::narrow<ioda::Dimensions_t>(numObjects))
        throw Exception("Size mismatch between underlying object and user-provided data range.",
                        ioda_Here());

      detail::PointerOwner pointerOwner = getTypeProvider()->getReturnedPointerOwner();
      Marshaller m(pointerOwner);
      auto p = m.prep_deserialize(numObjects);
      read(gsl::make_span<char>(reinterpret_cast<char*>(p->DataPointers.data()),
                                p->DataPointers.size() * Marshaller::bytesPerElement_),
           TypeWrapper::GetType(getTypeProvider()));
      m.deserialize(p, data);

      return Attribute_Implementation{backend_};
    } catch (...) {
      std::throw_with_nested(Exception(ioda_Here()));
    }
  }

  /// \brief Vector read convenience function.
  /// \tparam DataType is the type of the data. I.e. float, int, int32_t, uint16_t, std::string, etc.
  /// \param data is a vector acting as a data buffer that is filled with the metadata's contents.
  ///   It gets resized as needed.
  /// \returns Another instance of this Attribute. Used for operation chaining.
  /// \note data will be stored in row-major order.
  template <class DataType>
  Attribute_Implementation read(std::vector<DataType>& data) const {
    data.resize(getDimensions().numElements);
    return read<DataType>(gsl::make_span<DataType>(data.data(), data.size()));
  }

  /// \brief Valarray read convenience function.
  /// \tparam DataType is the type of the data. I.e. float, int, int32_t, uint16_t, std::string, etc.
  /// \param data is a valarray acting as a data buffer that is filled with the metadata's contents.
  ///   It gets resized as needed.
  /// \returns Another instance of this Attribute. Used for operation chaining.
  /// \note data will be stored in row-major order.
  template <class DataType>
  Attribute_Implementation read(std::valarray<DataType>& data) const {
    data.resize(getDimensions().numElements);
    return read<DataType>(gsl::make_span<DataType>(std::begin(data), std::end(data)));
  }

  /// \brief Read into a single value (convenience function).
  /// \tparam DataType is the type of the data. I.e. float, int, int32_t, uint16_t, std::string, etc.
  /// \param data is where the datum is read to.
  /// \returns Another instance of this Attribute. Used for operation chaining.
  /// \throws ioda::Exception if the underlying data have multiple elements.
  template <class DataType>
  Attribute_Implementation read(DataType& data) const {
    try {
      if (getDimensions().numElements != 1)
        throw Exception("Wrong number of elements. Use a different read() method.", ioda_Here());
      return read<DataType>(gsl::make_span<DataType>(&data, 1));
    } catch (...) {
      std::throw_with_nested(Exception(ioda_Here()));
    }
  }

  /// \brief Read a single value (convenience function).
  /// \tparam DataType is the type of the data. I.e. float, int, int32_t, uint16_t, std::string, etc.
  /// \returns A datum of type DataType.
  /// \throws ioda::Exception if the underlying data have size greater than 1.
  /// \note The Python function is read_datum_*
  template <class DataType>
  DataType read() const {
    DataType ret;
    read<DataType>(ret);
    return ret;
  }

  /// \brief Read into a new vector. Python convenience function.
  /// \tparam DataType is the type of the data.
  /// \note The Python function is read_list_*
  template <class DataType>
  std::vector<DataType> readAsVector() const {
    std::vector<DataType> data(getDimensions().numElements);
    read<DataType>(gsl::make_span<DataType>(data.data(), data.size()));
    return data;
  }

  /// \brief Read data into an Eigen::Array, Eigen::Matrix, Eigen::Map, etc.
  /// \tparam EigenClass is a template pointing to the Eigen object.
  ///   This template must provide the EigenClass::Scalar typedef.
  /// \tparam Resize indicates whether the Eigen object should be resized
  ///   if there is a dimension mismatch. Not all Eigen objects can be resized.
  /// \param res is the Eigen object.
  /// \returns Another instance of this Attribute. Used for operation chaining.
  /// \throws ioda::Exception if the attribute's dimensionality is
  ///   too high.
  /// \throws ioda::Exception if resize = false and there is a dimension mismatch.
  /// \note When reading in a 1-D object, the data are read as a column vector.
  template <class EigenClass, bool Resize = detail::EigenCompat::CanResize<EigenClass>::value>
  Attribute_Implementation readWithEigenRegular(EigenClass& res) const {
#if 1  //__has_include("Eigen/Dense")
    try {
      typedef typename EigenClass::Scalar ScalarType;

      static_assert(
        !(Resize && !detail::EigenCompat::CanResize<EigenClass>::value),
        "This object cannot be resized, but you have specified that a resize is required.");

      // Check that the dimensionality is 1 or 2.
      const auto dims = getDimensions();
      if (dims.dimensionality > 2)
        throw Exception(
          "Dimensionality too high for a regular Eigen read. Use "
          "Eigen::Tensor reads instead.",
          ioda_Here());

      int nDims[2] = {1, 1};
      if (dims.dimsCur.size() >= 1) nDims[0] = gsl::narrow<int>(dims.dimsCur[0]);
      if (dims.dimsCur.size() >= 2) nDims[1] = gsl::narrow<int>(dims.dimsCur[1]);

      // Resize if needed.
      if (Resize)
        detail::EigenCompat::DoEigenResize(res, nDims[0],
                                           nDims[1]);  // nullop if the size is already correct.
      else if (dims.numElements != (size_t)(res.rows() * res.cols()))
        throw Exception("Size mismatch", ioda_Here());

      // Array copy to preserve row vs column major format.
      // Should be optimized away by the compiler if unneeded.
      // Note to the reader: We are reading in the data to a temporary object.
      // We can size _this_ temporary object however we want.
      // The temporary is used to swap row / column indices if needed.
      // It should be optimized away if not needed... making sure this happens is a todo.
      Eigen::Array<ScalarType, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> data_in(res.rows(),
                                                                                res.cols());

      auto ret = read<ScalarType>(gsl::span<ScalarType>(data_in.data(), dims.numElements));
      res      = data_in;
      return ret;
    } catch (...) {
      std::throw_with_nested(Exception(ioda_Here()));
    }
#else
    static_assert(false, "The Eigen headers cannot be found, so this function cannot be used.");
#endif
  }

  /// \brief Read data into an Eigen::Array, Eigen::Matrix, Eigen::Map, etc.
  /// \tparam EigenClass is a template pointing to the Eigen object.
  ///   This template must provide the EigenClass::Scalar typedef.
  /// \param res is the Eigen object.
  /// \returns Another instance of this Attribute. Used for operation chaining.
  /// \throws ioda::Exception if there is a size mismatch.
  /// \note When reading in a 1-D object, the data are read as a column vector.
  template <class EigenClass>
  Attribute_Implementation readWithEigenTensor(EigenClass& res) const {
#if 1  //__has_include("unsupported/Eigen/CXX11/Tensor")
    try {
      // Check dimensionality of source and destination
      const auto ioda_dims  = getDimensions();
      const auto eigen_dims = ioda::detail::EigenCompat::getTensorDimensions(res);
      if (ioda_dims.numElements != eigen_dims.numElements)
        throw Exception("Size mismatch for Eigen Tensor-like read.", ioda_Here());

      auto sp = (gsl::make_span(res.data(), eigen_dims.numElements));
      return read(sp);
    } catch (...) {
      std::throw_with_nested(Exception(ioda_Here()));
    }
#else
    static_assert(
      false, "The Eigen unsupported/ headers cannot be found, so this function cannot be used.");
#endif
  }

  /// \internal Python binding function
  template <class EigenClass>
  EigenClass _readWithEigenRegular_python() const {
    EigenClass data;
    readWithEigenRegular(data);
    return data;
  }

  /// @}
  /// @name Type-querying Functions
  /// @{

  /// \brief Get Attribute type.
  virtual Type getType() const;
  /// \brief Get Attribute type.
  inline Type type() const { return getType(); }

  /// Query the backend and get the type provider.
  virtual detail::Type_Provider* getTypeProvider() const;

  /// \brief Convenience function to check an Attribute's storage type.
  /// \tparam DataType is the type of the data. I.e. float, int, int32_t, uint16_t, std::string, etc.
  /// \returns True if the type matches
  /// \returns False (0) if the type does not match
  /// \throws ioda::Exception if an error occurred.
  template <class DataType>
  bool isA() const {
    Type templateType = Types::GetType_Wrapper<DataType>::GetType(getTypeProvider());

    return isA(templateType);
  }
  /// Hand-off to the backend to check equivalence
  virtual bool isA(Type lhs) const;

  /// Python compatability function
  inline bool isA(BasicTypes dataType) { return isA(Type(dataType, getTypeProvider())); }
  /// \internal pybind11
  inline bool _py_isA2(BasicTypes dataType) { return isA(dataType); }

  /// @}
  /// @name Data Space-Querying Functions
  /// @{

  // Get an Attribute's dataspace
  // virtual Encapsulated_Handle getSpace() const;

  /// \brief Get Attribute's dimensions.
  virtual Dimensions getDimensions() const;

  /// @}
};

// extern template class Attribute_Base<>;
}  // namespace detail

/** \brief This class represents attributes, which may be attached to both Variables and Groups.
 *  \ingroup ioda_cxx_attribute
 *
 * Attributes are used to store small objects that get tagged to a Variable or a
 * Group to provide context to users and other programs. Attributes include
 * descriptions, units, alternate names, dimensions, and similar constructs.
 * Attributes may have different types (ints, floats, datetimes, strings, etc.),
 * and may be 0- or 1-dimensional.
 *
 * We can open an Attribute from a Has_Attribute object, which is a member of Groups and
 * Variables.
 *
 * \note Multidimensional attributes are supported by some of the underlying backends,
 * like HDF5, but are incompatible with the NetCDF file format.
 * \see Has_Attribute for the class that can create and open new Attribute objects.
 * \throws ioda::Exception on all exceptions.
 **/
class IODA_DL Attribute : public detail::Attribute_Base<> {
public:
  Attribute();
  Attribute(std::shared_ptr<detail::Attribute_Backend> b);
  Attribute(const Attribute&);
  Attribute& operator=(const Attribute&);
  virtual ~Attribute();

  /// @name Python compatability objects
  /// @{

  detail::python_bindings::AttributeIsA<Attribute> _py_isA;

  detail::python_bindings::AttributeReadSingle<Attribute> _py_readSingle;
  detail::python_bindings::AttributeReadVector<Attribute> _py_readVector;
  detail::python_bindings::AttributeReadNPArray<Attribute> _py_readNPArray;

  detail::python_bindings::AttributeWriteSingle<Attribute> _py_writeSingle;
  detail::python_bindings::AttributeWriteVector<Attribute> _py_writeVector;
  detail::python_bindings::AttributeWriteNPArray<Attribute> _py_writeNPArray;

  /// @}
};

namespace detail {
/// \brief Attribute backends inherit from this.
class IODA_DL Attribute_Backend : public Attribute_Base<> {
public:
  virtual ~Attribute_Backend();

protected:
  Attribute_Backend();
};
}  // namespace detail
}  // namespace ioda

/// @} // End Doxygen block

max_stars_repo_name: NOAA-EMC/ioda

max_issues_repo_name: NOAA-EMC/ioda

max_forks_repo_name: NOAA-EMC/ioda