// (c) 2008 Steven Gratton
// Guided in part by samples from the AMD CAL SDK
// and currently links with their Timer code

#include <iostream>
#include <iomanip>
#include <string>

#include <time.h>
#include <cmath>

#include "cal.h"
#include "calcl.h"

#include "bothstripedmmm.h"

#include "Timer.h"

#define PERFCOUNT

#ifdef PERFCOUNT
#include "cal_ext.h"
#include "cal_ext_counter.h"
#endif


std::string ILcheck=
  "il_ps_2_0\n"
  "ret_dyn\n"
  "end\n";

std::string ILcheck2=
  "il_ps_2_0\n"
  "dcl_input_position_interp(linear_noperspective) vWinCoord0.xy\n"
  "dcl_resource_id(0)_type(2d,unnorm)_fmtx(float)_fmty(float)_fmtz(float)_fmtw(float)\n"
  "dcl_resource_id(1)_type(2d,unnorm)_fmtx(float)_fmty(float)_fmtz(float)_fmtw(float)\n"
  "dcl_cb cb0[1] ; k4max, pitch of mat2 as floats \n"
  "sample_resource(0)_sampler(0) r0, vWinCoord0.xy\n"
  "sample_resource(1)_sampler(0) r1, vWinCoord0.xy\n"
  "imad r2.x___,vWinCoord0.y,cb0[0].y,vWinCoord0.x \n"
  "mov g[0],vWinCoord0.xy11\n"
  "ret_dyn\n"
  "end\n";

static CPerfCounter timer;

using namespace std;

#ifdef PERFCOUNT
static PFNCALCTXCREATECOUNTER  calCtxCreateCounterExt;
static PFNCALCTXDESTROYCOUNTER calCtxDestroyCounterExt;
static PFNCALCTXBEGINCOUNTER   calCtxBeginCounterExt;
static PFNCALCTXENDCOUNTER     calCtxEndCounterExt;
static PFNCALCTXGETCOUNTER     calCtxGetCounterExt;
#endif

void cpumatmult(int m,int k,int n,float* a,float* b,float* c)
{

  printf("*on the cpu...*\n");
  int i,j,p;
  float sum;
  volatile clock_t internalcpu;

  internalcpu=clock();
  for (i=0;i<m;i++)
    {
      for (j=0;j<n;j++)
        {
	  float tmp=0.f; 
	  for (p=0;p<k;p++)
	    {
	      tmp+=a[i*k+p]*b[p*n+j];
	    }
	  c[i*k+j]=tmp;
        }
    }
  internalcpu=clock()-internalcpu;
  printf("Internal cputime=%u.\n",internalcpu);
}

void maketestmat(int height,int width,float* mat)
{
  for (int i=0;i<height;i++){
    for (int j=0;j<width;j++){
      mat[i*width+j]=(float) (i*width+j);
    }
  }
}

void makediagtestmat(int height,int width,float* mat)
{
  for (int i=0;i<height;i++){
    for (int j=0;j<width;j++){
      mat[i*width+j]=0.f;
    }
    mat[i*width+i]=(float) i;
  }

}


void makeidentitymat(int height,int width,float* mat)
{
  for (int i=0;i<height;i++){
    for (int j=0;j<width;j++){
      mat[i*width+j]=0.f;
    }
    mat[i*width+i]=1.f;
  }
}


void  makestripe(int h,int w,float* mat,float* s0,float* s1,float* s2,float* s3,
		 float* s4,float* s5,float* s6,float* s7)
{
  for (int i=0;i<h;i++){
    int ii=i/8;
    int stripe=i%8;
    switch (stripe)
      {
      case 0:
    for (int j=0;j<w;j++){
      s0[ii*w+j]=mat[i*w+j];
    }
    break;
      case 1:
    for (int j=0;j<w;j++){
      s1[ii*w+j]=mat[i*w+j];
    }
    break;
      case 2:
    for (int j=0;j<w;j++){
      s2[ii*w+j]=mat[i*w+j];
    }
    break;
      case 3:
    for (int j=0;j<w;j++){
      s3[ii*w+j]=mat[i*w+j];
    }
    break;
      case 4:
    for (int j=0;j<w;j++){
      s4[ii*w+j]=mat[i*w+j];
    }
    break;
      case 5:
    for (int j=0;j<w;j++){
      s5[ii*w+j]=mat[i*w+j];
    }
    break;
      case 6:
    for (int j=0;j<w;j++){
      s6[ii*w+j]=mat[i*w+j];
    }
    break;
      case 7:
    for (int j=0;j<w;j++){
      s7[ii*w+j]=mat[i*w+j];
    }
    break;
      }
  }
}

void make1stripe(int h, int w, float* mat,
		 float* s00, float* s01, float* s02, float* s03,
		 float* s10, float* s11, float* s12, float* s13,
		 float* s20, float* s21, float* s22, float* s23,
		 float* s30, float* s31, float* s32, float* s33)
{
  int w4=w/4;
  int w16=w/16;
  for(int i=0; i<h; i++){
    int ii=i/4;
    int stripe=i%4;
    switch (stripe) {
    case 0:
      for (int j16=0;j16<w16;j16++){
	s00[ii*w4+4*j16]=mat[i*w+16*j16];
	s00[ii*w4+4*j16+1]=mat[i*w+16*j16+1];
	s00[ii*w4+4*j16+2]=mat[i*w+16*j16+2];
	s00[ii*w4+4*j16+3]=mat[i*w+16*j16+3];
	s01[ii*w4+4*j16]=mat[i*w+16*j16+4];
	s01[ii*w4+4*j16+1]=mat[i*w+16*j16+5];
	s01[ii*w4+4*j16+2]=mat[i*w+16*j16+6];
	s01[ii*w4+4*j16+3]=mat[i*w+16*j16+7];
	s02[ii*w4+4*j16]=mat[i*w+16*j16+8];
	s02[ii*w4+4*j16+1]=mat[i*w+16*j16+9];
	s02[ii*w4+4*j16+2]=mat[i*w+16*j16+10];
	s02[ii*w4+4*j16+3]=mat[i*w+16*j16+11];
	s03[ii*w4+4*j16]=mat[i*w+16*j16+12];
	s03[ii*w4+4*j16+1]=mat[i*w+16*j16+13];
	s03[ii*w4+4*j16+2]=mat[i*w+16*j16+14];
	s03[ii*w4+4*j16+3]=mat[i*w+16*j16+15];
      }
      break;
    case 1:
      for (int j16=0;j16<w16;j16++){
	s10[ii*w4+4*j16]=mat[i*w+16*j16];
	s10[ii*w4+4*j16+1]=mat[i*w+16*j16+1];
	s10[ii*w4+4*j16+2]=mat[i*w+16*j16+2];
	s10[ii*w4+4*j16+3]=mat[i*w+16*j16+3];
	s11[ii*w4+4*j16]=mat[i*w+16*j16+4];
	s11[ii*w4+4*j16+1]=mat[i*w+16*j16+5];
	s11[ii*w4+4*j16+2]=mat[i*w+16*j16+6];
	s11[ii*w4+4*j16+3]=mat[i*w+16*j16+7];
	s22[ii*w4+4*j16]=mat[i*w+16*j16+8];
	s22[ii*w4+4*j16+1]=mat[i*w+16*j16+9];
	s22[ii*w4+4*j16+2]=mat[i*w+16*j16+10];
	s22[ii*w4+4*j16+3]=mat[i*w+16*j16+11];
	s33[ii*w4+4*j16]=mat[i*w+16*j16+12];
	s33[ii*w4+4*j16+1]=mat[i*w+16*j16+13];
	s33[ii*w4+4*j16+2]=mat[i*w+16*j16+14];
	s33[ii*w4+4*j16+3]=mat[i*w+16*j16+15];
      }
      break;
    case 2:
      for (int j16=0;j16<w16;j16++){
	s20[ii*w4+4*j16]=mat[i*w+16*j16];
	s20[ii*w4+4*j16+1]=mat[i*w+16*j16+1];
	s20[ii*w4+4*j16+2]=mat[i*w+16*j16+2];
	s20[ii*w4+4*j16+3]=mat[i*w+16*j16+3];
	s21[ii*w4+4*j16]=mat[i*w+16*j16+4];
	s21[ii*w4+4*j16+1]=mat[i*w+16*j16+5];
	s21[ii*w4+4*j16+2]=mat[i*w+16*j16+6];
	s21[ii*w4+4*j16+3]=mat[i*w+16*j16+7];
	s22[ii*w4+4*j16]=mat[i*w+16*j16+8];
	s22[ii*w4+4*j16+1]=mat[i*w+16*j16+9];
	s22[ii*w4+4*j16+2]=mat[i*w+16*j16+10];
	s22[ii*w4+4*j16+3]=mat[i*w+16*j16+11];
	s23[ii*w4+4*j16]=mat[i*w+16*j16+12];
	s23[ii*w4+4*j16+1]=mat[i*w+16*j16+13];
	s23[ii*w4+4*j16+2]=mat[i*w+16*j16+14];
	s23[ii*w4+4*j16+3]=mat[i*w+16*j16+15];
      }
      break;
    case 3:
      for (int j16=0;j16<w16;j16++){
	s30[ii*w4+4*j16]=mat[i*w+16*j16];
	s30[ii*w4+4*j16+1]=mat[i*w+16*j16+1];
	s30[ii*w4+4*j16+2]=mat[i*w+16*j16+2];
	s30[ii*w4+4*j16+3]=mat[i*w+16*j16+3];
	s31[ii*w4+4*j16]=mat[i*w+16*j16+4];
	s31[ii*w4+4*j16+1]=mat[i*w+16*j16+5];
	s31[ii*w4+4*j16+2]=mat[i*w+16*j16+6];
	s31[ii*w4+4*j16+3]=mat[i*w+16*j16+7];
	s32[ii*w4+4*j16]=mat[i*w+16*j16+8];
	s32[ii*w4+4*j16+1]=mat[i*w+16*j16+9];
	s32[ii*w4+4*j16+2]=mat[i*w+16*j16+10];
	s32[ii*w4+4*j16+3]=mat[i*w+16*j16+11];
	s33[ii*w4+4*j16]=mat[i*w+16*j16+12];
	s33[ii*w4+4*j16+1]=mat[i*w+16*j16+13];
	s33[ii*w4+4*j16+2]=mat[i*w+16*j16+14];
	s33[ii*w4+4*j16+3]=mat[i*w+16*j16+15];
      }
      break;
    }
  }
}

void  destripe(int h,int w,float* mat,float* s0,float* s1,float* s2,float* s3,
		 float* s4,float* s5,float* s6,float* s7)
{
  for (int i=0;i<h;i++){
    int ii=i/8;
    int stripe=i%8;
    switch (stripe)
      {
      case 0:
    for (int j=0;j<w;j++){
      mat[i*w+j]=s0[ii*w+j];
    }
    break;
      case 1:
    for (int j=0;j<w;j++){
      mat[i*w+j]=s1[ii*w+j];
    }
    break;
      case 2:
    for (int j=0;j<w;j++){
      mat[i*w+j]=s2[ii*w+j];
    }
    break;
      case 3:
    for (int j=0;j<w;j++){
      mat[i*w+j]=s3[ii*w+j];
    }
    break;
      case 4:
    for (int j=0;j<w;j++){
      mat[i*w+j]=s4[ii*w+j];
    }
    break;
      case 5:
    for (int j=0;j<w;j++){
      mat[i*w+j]=s5[ii*w+j];
    }
    break;
      case 6:
    for (int j=0;j<w;j++){
      mat[i*w+j]=s6[ii*w+j];
    }
    break;
      case 7:
    for (int j=0;j<w;j++){
      mat[i*w+j]=s7[ii*w+j];
    }
    break;
      }
  }
}

void dispmat(int height,int width,float* mat)
{
  cout << "In dispmat." << endl;
  for (int i=0;i<height;i++){
    for (int j=0;j<width;j++){
      printf("%10.2e ",mat[i*width+j]);
      if (!((j+1)&0x3)) printf("  ");
    }
    printf("\n");
    if (!((i+1)&0x3)) printf("\n");
  }
  printf("\n");
}

void dispmatpart(int height,int width,int maximum,float* mat)
{
  int viewheight=height;
  if (maximum<viewheight) viewheight=maximum;
  int viewwidth=width;
  if (maximum<viewwidth) viewwidth=maximum;


  cout << "In dispmat." << endl;
  for (int i=0;i<viewheight;i++){
    for (int j=0;j<viewwidth;j++){
      printf("%12.4e ",mat[i*width+j]);
      if (!((j+1)&0x3)) printf("  ");
    }
    printf("\n");
    if (!((i+1)&0x3)) printf("\n");
  }
  printf("\n");
}

void dispdiag(int height,int width,float* mat)
{
  int top=0;
  if (width<=height) 
    top=width;
  else
    top=height;
  for (int i=0;i<top;i++){
    printf("%10.6e ",mat[i*width+i]);
    if (!((i+1)&0x3)) printf("  ");
    if (!((i+1)&0x1f)) printf("\n");
  }
  printf("\n");
}

void dispdiagpart(int height,int width,int maximum, float* mat)
{
  int top=height;
  if (width<height) top=width;
  if (maximum<top) top=maximum;
  for (int i=0;i<top;i++){
    printf("%10.6e ",mat[i*width+i]);
    if (!((i+1)&0x3)) printf("  ");
    if (!((i+1)&0x1f)) printf("\n");
  }
  printf("\n");
}

void dispmattr(int height,int width,float* mat)
{
  for (int i=0;i<height;i++){
    for (int j=0;j<width;j++){
      printf("%10.6f ",mat[j*width+i]);
      if (!((j+1)&0x3)) printf("  ");
    }
    printf("\n");
    if (!((i+1)&0x3)) printf("\n");
  }
  printf("\n");
}


int copytogpu(int width,int height,float* cpumat,CALresource gpumat)
{
  float* gpuptr=NULL;
  CALuint gpupitch=0;
  calResMap((CALvoid**)&gpuptr, &gpupitch, gpumat, 0);
  cout << "pitch=" << gpupitch << endl;
  for (int i = 0; i < height; ++i)
    {
      float* tmp = &gpuptr[i * gpupitch*4];
      for(int j = 0; j < width; ++j)
	{
	  tmp[j]=cpumat[i*width+j];
	}
    }
  calResUnmap(gpumat);
  return (int) gpupitch;
}

int zeroongpu(int width,int height,CALresource gpumat)
{
  float* gpuptr=NULL;
  CALuint gpupitch=0;
  calResMap((CALvoid**)&gpuptr, &gpupitch, gpumat, 0);
  memset(gpuptr,0,height*gpupitch*4*sizeof(float));
  calResUnmap(gpumat);
  return (int) gpupitch;
}

void copytocpu(int width,int height,float* cpumat,CALresource gpumat)
{
  cout << "In copytocpu." << endl;
  float* gpuptr=NULL;
  CALuint gpupitch=0;
  calResMap((CALvoid**)&gpuptr, &gpupitch, gpumat, 0);
  for (int i = 0; i < height; ++i)
    {
      float* tmp = &gpuptr[i * gpupitch*4];
      for(int j = 0; j < width; ++j)
	{
	  cpumat[i*width+j]=tmp[j];
	}
    }
  calResUnmap(gpumat);
}

int main(int argc, char** argv)
{
  //  multiply an m*k matrix by a k*n one to get a m*n one
  //  where m*k means m rows and k columns
  //  i.e. an array of width k and height m
  //  All arrays are stored row major
  //  It can get a bit confusing because "graphics-oriented"
  //  calls typically write width, height
  //  whereas matrix ones do height,width...

  int m=4096,k=4096,n=4096;
  int k4=k/4;
  int n4=n/4;
  int m4=m/4;
  int m8=m/8;
  int k8=k/8;
  int n8=n/8;
  int n16=n/16;

  //NEED TO CHECK k size...
  if((m%8)!=0||(k%16)!=0||(n%16)!=0) 
    {
      cout << "Sorry, matrices aren't appropriately sized." << endl; 
      return -1;
    };

  float* cpumat0=new float[m*k];
  float* cpumat1=new float[k*n];
  float* cpumat2=new float[m*n];

  float* cpumat0s0=new float[m8*k];
  float* cpumat0s1=new float[m8*k];
  float* cpumat0s2=new float[m8*k];
  float* cpumat0s3=new float[m8*k];
  float* cpumat0s4=new float[m8*k];
  float* cpumat0s5=new float[m8*k];
  float* cpumat0s6=new float[m8*k];
  float* cpumat0s7=new float[m8*k];

  float* cpumat1s00=new float[k4*n4];
  float* cpumat1s01=new float[k4*n4];
  float* cpumat1s02=new float[k4*n4];
  float* cpumat1s03=new float[k4*n4];
  float* cpumat1s10=new float[k4*n4];
  float* cpumat1s11=new float[k4*n4];
  float* cpumat1s12=new float[k4*n4];
  float* cpumat1s13=new float[k4*n4];
  float* cpumat1s20=new float[k4*n4];
  float* cpumat1s21=new float[k4*n4];
  float* cpumat1s22=new float[k4*n4];
  float* cpumat1s23=new float[k4*n4];
  float* cpumat1s30=new float[k4*n4];
  float* cpumat1s31=new float[k4*n4];
  float* cpumat1s32=new float[k4*n4];
  float* cpumat1s33=new float[k4*n4];

  // a "matrix" call so height, width...
  makeidentitymat(m,k,cpumat0);
  maketestmat(k,n,cpumat1);

dispmatpart(m,k,16,cpumat0);

  makestripe(m,k,cpumat0,cpumat0s0,cpumat0s1,cpumat0s2,cpumat0s3,cpumat0s4
	     ,cpumat0s5,cpumat0s6,cpumat0s7);

  dispmatpart(m8,k,16,cpumat0s0);
  dispmatpart(m8,k,16,cpumat0s1);

  
dispmatpart(k,n,16,cpumat1);
  make1stripe(k,n,cpumat1,
	      cpumat1s00,cpumat1s01,cpumat1s02,cpumat1s03,
	      cpumat1s10,cpumat1s11,cpumat1s12,cpumat1s13,
	      cpumat1s20,cpumat1s21,cpumat1s22,cpumat1s23,
	      cpumat1s30,cpumat1s31,cpumat1s32,cpumat1s33);

  dispmatpart(k4,n4,16,cpumat1s00);
  dispmatpart(k4,n4,16,cpumat1s01);
  dispmatpart(k4,n4,16,cpumat1s10);

  /*
    float* cpumat3=new float[m*n];
    cpumatmult(m,k,n,cpumat0,cpumat1,cpumat2);
    dispmat(m,n,cpumat2);
    delete[] cpumat3;
  */

  std::string kernel0 = fastmult;
  // kernel0=ILcheck2;

  calInit();
  CALuint numDevices = 0;
  calDeviceGetCount(&numDevices);

  cout << "Num devices =" << numDevices << endl;

  CALdevice device = 0;
  calDeviceOpen(&device, 0);

  CALdeviceinfo info;
  calDeviceGetInfo(&info, 0);


#ifdef PERFCOUNT
  if (calExtSupported((CALextid)CAL_EXT_COUNTERS) != CAL_RESULT_OK)
    {
      return 1;
    }
        
  if (calExtGetProc((CALextproc*)&calCtxCreateCounterExt, (CALextid)CAL_EXT_COUNTERS, "calCtxCreateCounter"))
    {
      return 1;
    }

  if (calExtGetProc((CALextproc*)&calCtxDestroyCounterExt, (CALextid)CAL_EXT_COUNTERS, "calCtxDestroyCounter"))
    {
      return 1;
    }
    
  if (calExtGetProc((CALextproc*)&calCtxBeginCounterExt, (CALextid)CAL_EXT_COUNTERS, "calCtxBeginCounter"))
    {
      return 1;
    }
    
  if (calExtGetProc((CALextproc*)&calCtxEndCounterExt, (CALextid)CAL_EXT_COUNTERS, "calCtxEndCounter"))
    {
      return 1;
    }

  if (calExtGetProc((CALextproc*)&calCtxGetCounterExt, (CALextid)CAL_EXT_COUNTERS, "calCtxGetCounter"))
    {
      return 1;
    }
#endif

  CALcontext ctx = 0;
  calCtxCreate(&ctx, device);

  CALobject obj0 = NULL;
  CALimage image0 = NULL;
  CALlanguage lang0 = CAL_LANGUAGE_IL;

  if (calclCompile(&obj0, lang0, kernel0.c_str(), info.target) != CAL_RESULT_OK)
    {
      fprintf(stdout, "Kernel0 compilation failed. Exiting.\n");
      cout << "Compile error: " << calclGetErrorString() << endl;
      return 1;
    }
  else
    {
      cout << "kernel0 compiled fine" << endl;
    };

  if (calclLink(&image0, &obj0, 1) != CAL_RESULT_OK)
    {
      fprintf(stdout, "Kernel0 linking failed. Exiting.\n");
      return 1;
    }

  // "graphics" so width, height...
  CALresource gpumat0s0=0;
  if(calResAllocLocal2D(&gpumat0s0, device,k4,m8, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat0 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat0 fine." << endl;
    }


  CALresource gpumat0s1=0;
  if(calResAllocLocal2D(&gpumat0s1, device,k4,m8, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat0 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat0 fine." << endl;
    }

  CALresource gpumat0s2=0;
  if(calResAllocLocal2D(&gpumat0s2, device,k4,m8, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat0 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat0 fine." << endl;
    }

  CALresource gpumat0s3=0;
  if(calResAllocLocal2D(&gpumat0s3, device,k4,m8, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat0 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat0 fine." << endl;
    }

  CALresource gpumat0s4=0;
  if(calResAllocLocal2D(&gpumat0s4, device,k4,m8, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat0 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat0 fine." << endl;
    }

  CALresource gpumat0s5=0;
  if(calResAllocLocal2D(&gpumat0s5, device,k4,m8, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat0 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat0 fine." << endl;
    }

  CALresource gpumat0s6=0;
  if(calResAllocLocal2D(&gpumat0s6, device,k4,m8, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat0 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat0 fine." << endl;
    }


  CALresource gpumat0s7=0;
  if(calResAllocLocal2D(&gpumat0s7, device,k4,m8, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat0 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat0 fine." << endl;
    }


  CALresource gpumat1s00=0;
  if(calResAllocLocal2D(&gpumat1s00, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s01=0;
  if(calResAllocLocal2D(&gpumat1s01, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }
  CALresource gpumat1s02=0;
  if(calResAllocLocal2D(&gpumat1s02, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s03=0;
  if(calResAllocLocal2D(&gpumat1s03, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s10=0;
  if(calResAllocLocal2D(&gpumat1s10, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s11=0;
  if(calResAllocLocal2D(&gpumat1s11, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s12=0;
  if(calResAllocLocal2D(&gpumat1s12, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s13=0;
  if(calResAllocLocal2D(&gpumat1s13, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s20=0;
  if(calResAllocLocal2D(&gpumat1s20, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s21=0;
  if(calResAllocLocal2D(&gpumat1s21, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s22=0;
  if(calResAllocLocal2D(&gpumat1s22, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s23=0;
  if(calResAllocLocal2D(&gpumat1s23, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s30=0;
  if(calResAllocLocal2D(&gpumat1s30, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }
  CALresource gpumat1s31=0;
  if(calResAllocLocal2D(&gpumat1s31, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s32=0;
  if(calResAllocLocal2D(&gpumat1s32, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }

  CALresource gpumat1s33=0;
  if(calResAllocLocal2D(&gpumat1s33, device, n16,k4, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("gpumat1 resource allocation failed.\n");
    }
  else
    {
      cout << "gpumat1 fine." << endl;
    }


  CALresource mat2=0;
  if(calResAllocLocal2D(&mat2, device, n4,m, 
			CAL_FORMAT_FLOAT_4, CAL_RESALLOC_GLOBAL_BUFFER)
     !=CAL_RESULT_OK) 
    {
      printf("mat2 resource allocation failed.\n");
    }
  else
    {
      cout << "mat2 fine." << endl;
    }


  CALresource const0=0;
  if(calResAllocLocal1D(&const0, device, 1, 
			CAL_FORMAT_FLOAT_4, 0)
     !=CAL_RESULT_OK) 
    {
      printf("const0 resource allocation failed.\n");
    }
  else
    {
      cout << "const0 fine." << endl;
    }

  // "graphics" so width, height...
  copytogpu(k,m8,cpumat0s0,gpumat0s0);
  copytogpu(k,m8,cpumat0s1,gpumat0s1);
  copytogpu(k,m8,cpumat0s2,gpumat0s2);
  copytogpu(k,m8,cpumat0s3,gpumat0s3);
  copytogpu(k,m8,cpumat0s4,gpumat0s4);
  copytogpu(k,m8,cpumat0s5,gpumat0s5);
  copytogpu(k,m8,cpumat0s6,gpumat0s6);
  copytogpu(k,m8,cpumat0s7,gpumat0s7);
  copytogpu(n16,k4,cpumat1s00,gpumat1s00);
  copytogpu(n16,k4,cpumat1s01,gpumat1s01);
  copytogpu(n16,k4,cpumat1s02,gpumat1s02);
  copytogpu(n16,k4,cpumat1s03,gpumat1s03);
  copytogpu(n16,k4,cpumat1s10,gpumat1s10);
  copytogpu(n16,k4,cpumat1s11,gpumat1s11);
  copytogpu(n16,k4,cpumat1s12,gpumat1s12);
  copytogpu(n16,k4,cpumat1s13,gpumat1s13);
  copytogpu(n16,k4,cpumat1s20,gpumat1s20);
  copytogpu(n16,k4,cpumat1s21,gpumat1s21);
  copytogpu(n16,k4,cpumat1s22,gpumat1s22);
  copytogpu(n16,k4,cpumat1s23,gpumat1s23);
  copytogpu(n16,k4,cpumat1s30,gpumat1s30);
  copytogpu(n16,k4,cpumat1s31,gpumat1s31);
  copytogpu(n16,k4,cpumat1s32,gpumat1s32);
  copytogpu(n16,k4,cpumat1s33,gpumat1s33);
  int pitch2=0;
  pitch2=zeroongpu(n,m,mat2);
  cout <<"pitch 2="<<pitch2<<"."<<endl;

  union {int i; float f;} tmpunion;
  tmpunion.i=pitch2;

  float* constdata=NULL;
  CALuint constpitch=0;     
  calResMap((void**)&constdata,&constpitch,const0,0);
  constdata[0]=(float) k4;
  constdata[1]=tmpunion.f;
  // dirty trick duplicating cb0.y and cb0.z in order to make compiler not keep
  // glob mem pointer registers kicking around from one side of the loop to the other...
  constdata[2]=tmpunion.f;

  calResUnmap(const0);

  CALmem mat0s0mem=0;
  CALmem mat0s1mem=0;
  CALmem mat0s2mem=0;
  CALmem mat0s3mem=0;
  CALmem mat0s4mem=0;
  CALmem mat0s5mem=0;
  CALmem mat0s6mem=0;
  CALmem mat0s7mem=0;
  CALmem mat1s00mem=0;
  CALmem mat1s01mem=0;
  CALmem mat1s02mem=0;
  CALmem mat1s03mem=0;
  CALmem mat1s10mem=0;
  CALmem mat1s11mem=0;
  CALmem mat1s12mem=0;
  CALmem mat1s13mem=0;
  CALmem mat1s20mem=0;
  CALmem mat1s21mem=0;
  CALmem mat1s22mem=0;
  CALmem mat1s23mem=0;
  CALmem mat1s30mem=0;
  CALmem mat1s31mem=0;
  CALmem mat1s32mem=0;
  CALmem mat1s33mem=0;
  CALmem mat2mem=0;
  CALmem const0mem=0;

  calCtxGetMem(&mat0s0mem, ctx, gpumat0s0);
  calCtxGetMem(&mat0s1mem, ctx, gpumat0s1);
  calCtxGetMem(&mat0s2mem, ctx, gpumat0s2);
  calCtxGetMem(&mat0s3mem, ctx, gpumat0s3);
  calCtxGetMem(&mat0s4mem, ctx, gpumat0s4);
  calCtxGetMem(&mat0s5mem, ctx, gpumat0s5);
  calCtxGetMem(&mat0s6mem, ctx, gpumat0s6);
  calCtxGetMem(&mat0s7mem, ctx, gpumat0s7);
  calCtxGetMem(&mat1s00mem, ctx, gpumat1s00);
  calCtxGetMem(&mat1s01mem, ctx, gpumat1s01);
  calCtxGetMem(&mat1s02mem, ctx, gpumat1s02);
  calCtxGetMem(&mat1s03mem, ctx, gpumat1s03);
  calCtxGetMem(&mat1s10mem, ctx, gpumat1s10);
  calCtxGetMem(&mat1s11mem, ctx, gpumat1s11);
  calCtxGetMem(&mat1s12mem, ctx, gpumat1s12);
  calCtxGetMem(&mat1s13mem, ctx, gpumat1s13);
  calCtxGetMem(&mat1s20mem, ctx, gpumat1s20);
  calCtxGetMem(&mat1s21mem, ctx, gpumat1s21);
  calCtxGetMem(&mat1s22mem, ctx, gpumat1s22);
  calCtxGetMem(&mat1s23mem, ctx, gpumat1s23);
  calCtxGetMem(&mat1s30mem, ctx, gpumat1s30);
  calCtxGetMem(&mat1s31mem, ctx, gpumat1s31);
  calCtxGetMem(&mat1s32mem, ctx, gpumat1s32);
  calCtxGetMem(&mat1s33mem, ctx, gpumat1s33);
  calCtxGetMem(&mat2mem, ctx, mat2);
  calCtxGetMem(&const0mem,ctx,const0);

  cout << "After GetMem, " << calGetErrorString() << "." << endl;

  CALmodule module0 = 0;

  calModuleLoad(&module0, ctx, image0);

  cout << "After ModuleLoad, " << calGetErrorString() << "." << endl;

  CALfunc func0 = 0;

  CALname matname0s0=0;
  CALname matname0s1=0;
  CALname matname0s2=0;
  CALname matname0s3=0;
  CALname matname0s4=0;
  CALname matname0s5=0;
  CALname matname0s6=0;
  CALname matname0s7=0;
  CALname matname1s00=0;
  CALname matname1s01=0;
  CALname matname1s02=0;
  CALname matname1s03=0;
  CALname matname1s10=0;
  CALname matname1s11=0;
  CALname matname1s12=0;
  CALname matname1s13=0;
  CALname matname1s20=0;
  CALname matname1s21=0;
  CALname matname1s22=0;
  CALname matname1s23=0;
  CALname matname1s30=0;
  CALname matname1s31=0;
  CALname matname1s32=0;
  CALname matname1s33=0;
  CALname matname2=0;
  CALname constname0=0;

  calModuleGetEntry(&func0, ctx, module0, "main");

  cout << "After GetEntry, " << calGetErrorString() << "." << endl;


  calModuleGetName(&matname0s0, ctx, module0, "i0");
  calModuleGetName(&matname0s1, ctx, module0, "i1");
  calModuleGetName(&matname0s2, ctx, module0, "i2");
  calModuleGetName(&matname0s3, ctx, module0, "i3");
  calModuleGetName(&matname0s4, ctx, module0, "i4");
  calModuleGetName(&matname0s5, ctx, module0, "i5");
  calModuleGetName(&matname0s6, ctx, module0, "i6");
  calModuleGetName(&matname0s7, ctx, module0, "i7");
  calModuleGetName(&matname1s00, ctx, module0, "i8");
  calModuleGetName(&matname1s01, ctx, module0, "i9");
  calModuleGetName(&matname1s02, ctx, module0, "i10");
  calModuleGetName(&matname1s03, ctx, module0, "i11");
  calModuleGetName(&matname1s10, ctx, module0, "i12");
  calModuleGetName(&matname1s11, ctx, module0, "i13");
  calModuleGetName(&matname1s12, ctx, module0, "i14");
  calModuleGetName(&matname1s13, ctx, module0, "i15");
  calModuleGetName(&matname1s20, ctx, module0, "i16");
  calModuleGetName(&matname1s21, ctx, module0, "i17");
  calModuleGetName(&matname1s22, ctx, module0, "i18");
  calModuleGetName(&matname1s23, ctx, module0, "i19");
  calModuleGetName(&matname1s30, ctx, module0, "i20");
  calModuleGetName(&matname1s31, ctx, module0, "i21");
  calModuleGetName(&matname1s32, ctx, module0, "i22");
  calModuleGetName(&matname1s33, ctx, module0, "i23");
  calModuleGetName(&matname2, ctx, module0, "g[]");
  calModuleGetName(&constname0, ctx, module0, "cb0");

  cout << "After GetName, " << calGetErrorString() << "." << endl;

  calCtxSetMem(ctx, matname0s0, mat0s0mem);
  calCtxSetMem(ctx, matname0s1, mat0s1mem);
  calCtxSetMem(ctx, matname0s2, mat0s2mem);
  calCtxSetMem(ctx, matname0s3, mat0s3mem);
  calCtxSetMem(ctx, matname0s4, mat0s4mem);
  calCtxSetMem(ctx, matname0s5, mat0s5mem);
  calCtxSetMem(ctx, matname0s6, mat0s6mem);
  calCtxSetMem(ctx, matname0s7, mat0s7mem);
  calCtxSetMem(ctx, matname1s00, mat1s00mem);
  calCtxSetMem(ctx, matname1s01, mat1s01mem);
  calCtxSetMem(ctx, matname1s02, mat1s02mem);
  calCtxSetMem(ctx, matname1s03, mat1s03mem);
  calCtxSetMem(ctx, matname1s10, mat1s10mem);
  calCtxSetMem(ctx, matname1s11, mat1s11mem);
  calCtxSetMem(ctx, matname1s12, mat1s12mem);
  calCtxSetMem(ctx, matname1s13, mat1s13mem);
  calCtxSetMem(ctx, matname1s20, mat1s20mem);
  calCtxSetMem(ctx, matname1s21, mat1s21mem);
  calCtxSetMem(ctx, matname1s22, mat1s22mem);
  calCtxSetMem(ctx, matname1s23, mat1s23mem);
  calCtxSetMem(ctx, matname1s30, mat1s30mem);
  calCtxSetMem(ctx, matname1s31, mat1s31mem);
  calCtxSetMem(ctx, matname1s32, mat1s32mem);
  calCtxSetMem(ctx, matname1s33, mat1s33mem);
  calCtxSetMem(ctx, matname2, mat2mem);
  calCtxSetMem(ctx, constname0, const0mem);

  cout << "After SetMem, " << calGetErrorString() << "." << endl;

  
  CALevent e = 0;
  CALevent f=0;

  timer.Reset();

  cout << "Just before running, " << calGetErrorString() << "." << endl;

  volatile clock_t gputime;

#ifdef PERFCOUNT

  CALcounter idleCounter;
  if (calCtxCreateCounterExt(&idleCounter, ctx, CAL_COUNTER_IDLE) != CAL_RESULT_OK)
    {
      return 1;
    }

  CALcounter cacheCounter;
  if (calCtxCreateCounterExt(&cacheCounter, ctx, CAL_COUNTER_INPUT_CACHE_HIT_RATE) != CAL_RESULT_OK)
    {
      return 1;
    }
    
  if (calCtxBeginCounterExt(ctx, idleCounter) != CAL_RESULT_OK)
    {
      return 1;
    }
    
  if (calCtxBeginCounterExt(ctx, cacheCounter) != CAL_RESULT_OK)
    {
      return 1;
    }
#endif

  // "graphics" so width, height...
  //CALdomain domain0 = {2,2,4,4};
  CALdomain domain0 = {0,0,n16,m8};
  //  CALdomain domain0 = {0,1,1,1};


  // xblock=16, yblock=64 seem best when looping over y first
  // xblock=32, yblock=16 seem best when looping over x first; seems best overall
       

  CALuint xstart=0;
  CALuint ystart=0;
  CALuint xblock=n16;
  CALuint yblock=m8;


  if ((n16%xblock)!=0||(m8%yblock)!=0) {cout << "block sizes no good." << endl; return 1;}

  //  testing from here...
  //  CALdomain domaintest={1,0,2,4};

  //calCtxFlush(ctx);

  //  calCtxRunProgram(&e, ctx, func0, &domaintest);
  //  while (calCtxIsEventDone(ctx, e) == CAL_RESULT_PENDING);e=0;


  //   cout << "After test, " << calGetErrorString() << endl;

  //  while (calCtxIsEventDone(ctx, e) == CAL_RESULT_PENDING);e=0;

  //  ... to here 

  gputime=clock();

  timer.Start();

  system("date");

  for (int kk=0;kk<50;kk++)
    {
      cout << clock()-gputime<< endl;
      cout << "kk=" << kk << endl;
      while(xstart<n16)
        {
          ystart=0;
          while(ystart<m8)
            {
              CALdomain domain1 = {xstart, ystart, (xstart+xblock), (ystart+yblock)};
	      //      if(((ystart+yblock)==m8)&&((xstart+xblock)==n16))
	      //	{	   
	      //calCtxRunProgram(&f, ctx, func0, &domain1);
	      //cout << "In f event." << endl;
	      //   calCtxFlush(ctx);

	      //	}
	      //else{
	      calCtxRunProgram(&e, ctx, func0, &domain1);
	      // calCtxFlush(ctx);
	      while (calCtxIsEventDone(ctx, e) == CAL_RESULT_PENDING);e=0;
	      //}

	      ystart+=yblock;
	      //	   		  cout << domain1.x << "," << domain1.y << "," <<  domain1.width << ","<< domain1.height << endl;
            };
	  xstart+=xblock;
	  //	  cout<< "ystart=" << ystart << endl;
        };

      //	      while (calCtxIsEventDone(ctx, f) == CAL_RESULT_PENDING);f=0;


      //  cout << "After spin, " << calGetErrorString() << "." << endl;


  
      //        calCtxRunProgram(&e, ctx, func0, &domain0);


      // cout << "At A, " << calGetErrorString() << "." << endl;

      //          while (calCtxIsEventDone(ctx, e) == CAL_RESULT_PENDING);e=0;
      //cout << "At B, " << calGetErrorString() << "." << endl;

      //system("date");

      xstart=0;
      ystart=0;

    };

  timer.Stop();

  gputime=clock()-gputime;

  cout << "After kernel0, " << calGetErrorString() << endl;

#ifdef PERFCOUNT
  if (calCtxEndCounterExt(ctx, idleCounter) != CAL_RESULT_OK)
    {
      return 1;
    }

  if (calCtxEndCounterExt(ctx, cacheCounter) != CAL_RESULT_OK)
    {
      return 1;
    }


  CALfloat idlePercentage = 0.0f;
  if (calCtxGetCounterExt(&idlePercentage, ctx, idleCounter) != CAL_RESULT_OK)
    {
      return 1;
    }

  CALfloat cachePercentage = 0.0f;
  if (calCtxGetCounterExt(&cachePercentage, ctx, cacheCounter) != CAL_RESULT_OK)
    {
      return 1;
    }
  printf("Idle percentage: %0.2f%% Cache hit rate: %0.2f%%\n", idlePercentage * 100.0f, cachePercentage * 100.0f);

  if (calCtxDestroyCounterExt(ctx, idleCounter) != CAL_RESULT_OK)
    {
      return 1;
    }

  if (calCtxDestroyCounterExt(ctx, cacheCounter) != CAL_RESULT_OK)
    {
      return 1;
    }

#endif

  cout << "gpu time=" << gputime << "/" << CLOCKS_PER_SEC << "=" <<gputime/1.e6f <<" s." <<endl;

  cout << "timer time=" << timer.GetElapsedTime() <<"s." << endl;

  system("date");

  cout << "After calculation, " << calGetErrorString() << "." << endl;

  // ambiguous, so chosen "graphics"...
  copytocpu(n,m,cpumat2,mat2);

  cout << "After copy, " << calGetErrorString() << "." << endl;

  // "matrix" so height, width...
  //dispmat(m,n,cpumat2);
  //dispdiag(m,n,cpumat2);
  // dispdiagpart(m,n,16,cpumat2);
  dispmatpart(m,n,64,cpumat2);

  cout << "After dispmat." << endl;

  delete[] cpumat0;
  delete[] cpumat0s0;
  delete[] cpumat0s1;
  delete[] cpumat0s2;
  delete[] cpumat0s3;
  delete[] cpumat0s4;
  delete[] cpumat0s5;
  delete[] cpumat0s6;
  delete[] cpumat0s7;
  delete[] cpumat1;
  delete[] cpumat1s00;
  delete[] cpumat1s01;
  delete[] cpumat1s02;
  delete[] cpumat1s03;
  delete[] cpumat1s10;
  delete[] cpumat1s11;
  delete[] cpumat1s12;
  delete[] cpumat1s13;
  delete[] cpumat1s20;
  delete[] cpumat1s21;
  delete[] cpumat1s22;
  delete[] cpumat1s23;
  delete[] cpumat1s30;
  delete[] cpumat1s31;
  delete[] cpumat1s32;
  delete[] cpumat1s33;
  delete[] cpumat2;
  
  calModuleUnload(ctx, module0);

  calclFreeImage(image0);
  calclFreeObject(obj0);

  calCtxReleaseMem(ctx,mat0s0mem);
  calCtxReleaseMem(ctx,mat0s1mem);
  calCtxReleaseMem(ctx,mat0s2mem);
  calCtxReleaseMem(ctx,mat0s3mem);
  calCtxReleaseMem(ctx,mat0s4mem);
  calCtxReleaseMem(ctx,mat0s5mem);
  calCtxReleaseMem(ctx,mat0s6mem);
  calCtxReleaseMem(ctx,mat0s7mem);
  calCtxReleaseMem(ctx,mat1s00mem);
  calCtxReleaseMem(ctx,mat1s01mem);
  calCtxReleaseMem(ctx,mat1s02mem);
  calCtxReleaseMem(ctx,mat1s03mem);
  calCtxReleaseMem(ctx,mat1s10mem);
  calCtxReleaseMem(ctx,mat1s11mem);
  calCtxReleaseMem(ctx,mat1s12mem);
  calCtxReleaseMem(ctx,mat1s13mem);
  calCtxReleaseMem(ctx,mat1s20mem);
  calCtxReleaseMem(ctx,mat1s21mem);
  calCtxReleaseMem(ctx,mat1s22mem);
  calCtxReleaseMem(ctx,mat1s23mem);
  calCtxReleaseMem(ctx,mat1s30mem);
  calCtxReleaseMem(ctx,mat1s31mem);
  calCtxReleaseMem(ctx,mat1s32mem);
  calCtxReleaseMem(ctx,mat1s33mem);
  calCtxReleaseMem(ctx,mat2mem);
  calCtxReleaseMem(ctx,const0mem);

  calResFree(gpumat0s0);
  calResFree(gpumat0s1);
  calResFree(gpumat0s2);
  calResFree(gpumat0s3);
  calResFree(gpumat0s4);
  calResFree(gpumat0s5);
  calResFree(gpumat0s6);
  calResFree(gpumat0s7);
  calResFree(gpumat1s00);
  calResFree(gpumat1s01);
  calResFree(gpumat1s02);
  calResFree(gpumat1s03);
  calResFree(gpumat1s10);
  calResFree(gpumat1s11);
  calResFree(gpumat1s12);
  calResFree(gpumat1s13);
  calResFree(gpumat1s20);
  calResFree(gpumat1s21);
  calResFree(gpumat1s22);
  calResFree(gpumat1s23);
  calResFree(gpumat1s30);
  calResFree(gpumat1s31);
  calResFree(gpumat1s32);
  calResFree(gpumat1s33);
  calResFree(mat2);
  calResFree(const0);

  calCtxDestroy(ctx);

  calDeviceClose(device);

  calShutdown();

  return 0;
}