vsldquantilerow.c

使用INTEL矢量统计类库的程序,包括以下功能: &#61623 Raw and central moments up to 4th order &#61623 Kurtosis and

/*******************************************************************************
!                             INTEL CONFIDENTIAL
!  Copyright(C) 2007-2008 Intel Corporation. All Rights Reserved.
!  The source code contained  or  described herein and all documents related to
!  the source code ("Material") are owned by Intel Corporation or its suppliers
!  or licensors.  Title to the  Material remains with  Intel Corporation or its
!  suppliers and licensors. The Material contains trade secrets and proprietary
!  and  confidential  information of  Intel or its suppliers and licensors. The
!  Material  is  protected  by  worldwide  copyright  and trade secret laws and
!  treaty  provisions. No part of the Material may be used, copied, reproduced,
!  modified, published, uploaded, posted, transmitted, distributed or disclosed
!  in any way without Intel's prior express written permission.
!  No license  under any  patent, copyright, trade secret or other intellectual
!  property right is granted to or conferred upon you by disclosure or delivery
!  of the Materials,  either expressly, by implication, inducement, estoppel or
!  otherwise.  Any  license  under  such  intellectual property  rights must be
!  express and approved by Intel in writing.
!
!*******************************************************************************
!  Content:
!    Calculation of quantiles and order statistics  Example Program Text
!******************************************************************************/


#include "mkl.h"
#include "vsl_ss.h"
#include <stdio.h>

#define RETURN_ON_ERROR                 \
    if(errcode<0)                        \
    {                                   \
        printf("Error: %i\n", errcode);  \
        printf("\nTEST FAILED\n");		\
        return errcode;                  \
    }

#define SEED    777
#define BRNG    VSL_BRNG_MCG31
#define METHOD  0

#define DIM 3      /* dimension of the task */ 
#define N   1000   /* number of observations */
#define M   11

int main()
{
    int i, j, errcode;
    int dim=DIM, n=N;
    int storage = VSL_SS_MATRIX_ROWS_STORAGE;
    int order_statistics_storage = VSL_SS_MATRIX_ROWS_STORAGE;
    VSLSSTaskPtr task=0;
    double xT[N][DIM];  /* matrix of observations */
    double x[DIM][N];  /* matrix of observations */
    double order_statistics[DIM][N]; /* matrix to store order statistics */
    double quant_order[M], quantiles[DIM][M];
    int quant_order_n=M;
    VSLStreamStatePtr stream;
    double a=0.0,sigma=1.0;
    int fail=0, more=0, less=0;

    for ( i = 0; i < M; i++ ) quant_order[i] = (double)i / (double)(M-1);

    /***** Initialize *****/
    errcode = vslNewStream( &stream, BRNG,  SEED );

    /***** Call RNG *****/
    errcode = vdRngGaussian( METHOD, stream, N*DIM, (double*)xT, a, sigma );

    /***** Deinitialize *****/
    errcode = vslDeleteStream( &stream );

    for ( j = 0; j < DIM; j++ )
    { 
        for ( i = 0; i < N; i++ ) 
        {
            x[j][i]=xT[i][j];
        }
    }
    

    /* Computation of deciles */  

    errcode = vsldSSNewTask( &task, &dim, &n, &storage, (double*)x, 0, 0 );
    RETURN_ON_ERROR;

    errcode = vsldSSEditQuantiles( task, &quant_order_n, quant_order,
                        (double*)quantiles, (double*)order_statistics, 
                        &order_statistics_storage );
    RETURN_ON_ERROR;

    errcode = vsldSSCompute(task, VSL_SS_QUANTILES|VSL_SS_ORDER_STATISTICS, 
                                 VSL_SS_FAST_METHOD );
    RETURN_ON_ERROR;

    errcode = vslSSDeleteTask( &task );
    RETURN_ON_ERROR;

    printf("1st 4 and last 4 observations in source matrix\n");
    for ( j = 0; j < DIM; j++ ) 
    {
        for ( i = 0; i < 4; i++ )
        {
            printf("%+.3f ", x[j][i] );
        }
        printf("     ...      ");
        for ( i = N-5; i < N; i++ )
        {
            printf("%+.3f ", x[j][i] );
        }
        printf("\n");
    }

    printf("Deciles of the observations for all variables:\n   ");
    for ( i = 0; i < M; i++ ) printf("D%i     ", i );
    printf("\n");
    for ( j = 0; j < DIM; j++ ) 
    {
        for ( i = 0; i < M; i++ ) printf("%+.3f ", quantiles[j][i] );
        printf("\n");
    }

    printf("1st 4 and last 4 observations in order statistics matrix\n");
    for ( j = 0; j < DIM; j++ ) 
    {
        for ( i = 0; i < 4; i++ )
        {
            printf("%+.3f ", order_statistics[j][i] );
        }
        printf("     ...      ");
        for ( i = N-5; i < N; i++ )
        {
            printf("%+.3f ", order_statistics[j][i] );
        }
        printf("\n");
    }

    for( j = 0; j < DIM; j++ )
    {
        for ( i = 0; i < N-1; i++ )
        {
            if(order_statistics[j][i] > order_statistics[j][i+1]) fail++;
        }
        
        more=less=0;

        for ( i = 0; i < N; i++ )
        {
            if(x[j][i] >= quantiles[j][M/2]) more++;
            if(x[j][i] <= quantiles[j][M/2]) less++;
        }

        if(more<(N/2+(N & 1)) || less<(N/2+(N & 1))) fail++;
    }

    printf("\nTEST PASSED\n");

    return 0;
}