C++ array2d类代码示例

int HMatrix::ParaDividing( 
		const array2d &vec2,
		array2d &product
	) const
{
	if(nColumns != (int)vec2[0].size())
	{
//		cout << "MatrixParaProduct has different columns for vec1 and vec2!" << endl;
		return -1;
	}

	product.resize( vec2.size(), array1d(nRows) );

	for (int i=0; i<(int)vec2.size(); i++)
	{
		for (int j=0; j<nRows; j++)
		{
			product[i][j] = 0.0;
			for (int k=0; k<nColumns; k++)
			{
				product[i][j] += data[j][k] / vec2[i][k];
			}
			if(product[i][j] == 0.0) product[i][j] = EPS;
		}
	}

	return 1;
}

    reference operator=(const array2d<float> arr){
      if(arr.width()!=width() || arr.height()!=height())
        throw "Error! Grids were not of the same size!";

      #pragma omp parallel for collapse(2)
      for(int x=0;x<width();x++)
      for(int y=0;y<height();y++)
        operator()(x,y)=arr(x,y);
    }

inline double
total(const array2d& u)
{
  double s = 0;
  const int nx = u.size_x();
  const int ny = u.size_y();
  for (int y = 1; y < ny - 1; y++)
    for (int x = 1; x < nx - 1; x++)
      s += u(x, y);
  return s;
}

inline void
time_step_indexed(array2d& u, const Constants& c)
{
  // compute du/dt
  array2d du(c.nx, c.ny, u.rate(), 0, u.cache_size());
  for (int y = 1; y < c.ny - 1; y++) {
    for (int x = 1; x < c.nx - 1; x++) {
      double uxx = (u(x - 1, y) - 2 * u(x, y) + u(x + 1, y)) / (c.dx * c.dx);
      double uyy = (u(x, y - 1) - 2 * u(x, y) + u(x, y + 1)) / (c.dy * c.dy);
      du(x, y) = c.dt * c.k * (uxx + uyy);
    }
  }
  // take forward Euler step
  for (uint i = 0; i < u.size(); i++)
    u[i] += du[i];
}

void d8_flow_directions(
  const array2d<T> &elevations,
  array2d<U> &flowdirs
){
  ProgressBar progress;

  diagnostic("Setting up the flow directions matrix...");
  flowdirs.copyprops(elevations);
  flowdirs.init(NO_FLOW);
  flowdirs.no_data=d8_NO_DATA;
  diagnostic("succeeded.\n");

  diagnostic("%%Calculating D8 flow directions...\n");
  progress.start( elevations.width()*elevations.height() );
  #pragma omp parallel for
  for(int x=0;x<elevations.width();x++){
    progress.update( x*elevations.height() );
    for(int y=0;y<elevations.height();y++)
      if(elevations(x,y)==elevations.no_data)
        flowdirs(x,y)=flowdirs.no_data;
      else
        flowdirs(x,y)=d8_FlowDir(elevations,x,y);
  }
  diagnostic_arg(SUCCEEDED_IN,progress.stop());
}

inline bool SaveGrayscaleToImageFile( const array2d<U8>& texel, const std::string& filepath )
{
	int x,y;
	int width  = texel.size_x();
	int height = texel.size_y();
	const int depth = 24;

	BitmapImage img( width, height, depth );

	for( y=0; y<height ; y++ )
	{
		for( x=0; x<width; x++ )
		{
			img.SetGrayscalePixel( x, y, texel(x,y) );
		}
	}

	return img.SaveToFile( filepath );
}

inline void
time_step_iterated(array2d& u, const Constants& c)
{
  // compute du/dt
  array2d du(c.nx, c.ny, u.rate(), 0, u.cache_size());
  for (typename array2d::iterator p = du.begin(); p != du.end(); p++) {
    int x = p.i();
    int y = p.j();
    if (1 <= x && x <= c.nx - 2 &&
        1 <= y && y <= c.ny - 2) {
      double uxx = (u(x - 1, y) - 2 * u(x, y) + u(x + 1, y)) / (c.dx * c.dx);
      double uyy = (u(x, y - 1) - 2 * u(x, y) + u(x, y + 1)) / (c.dy * c.dy);
      *p = c.dt * c.k * (uxx + uyy);
    }
  }
  // take forward Euler step
  for (typename array2d::iterator p = u.begin(), q = du.begin(); p != u.end(); p++, q++)
    *p += *q;
}

inline void GetPerlinTexture( const PerlinNoiseParams& params, array2d<float>& dest )
{
	Perlin pn( params.octaves, params.freq, params.amp, params.seed );

	float (*pPerlinFunc) (Perlin&,float,float,float,float);
	pPerlinFunc = params.tilable ? TilablePerlin : StdPerlin;

	const int w = dest.size_x();
	const int h = dest.size_y();
	float min_val = FLT_MAX, max_val = -FLT_MAX;
	for( int y=0; y<h; y++ )
	{
		for( int x=0; x<w; x++ )
		{
			float fx = (float)x / (float)w;
			float fy = (float)y / (float)h;
//			float f = pn.Get( fx, fy );
			float f = pPerlinFunc( pn, fx, fy, 1.0f, 1.0f );
			min_val = take_min( min_val, f );
			max_val = take_max( max_val, f );

			dest(x,y) = f;
		}
	}

	float val_range = max_val - min_val;

//	printf( "(min,max) = (%f,%f)\n", min_val, max_val );

	float dest_min   = params.min_value;
	float dest_range = params.max_value - params.min_value;
	for( int y=0; y<h; y++ )
	{
		for( int x=0; x<w; x++ )
		{
			float normalized_val = ( dest(x,y) - min_val ) / val_range; // normalized_val is [0,1]
			dest(x,y) = dest_min + normalized_val * dest_range;
		}
	}
}

int VMatrix::Dividing(
		const array2d &scaler,
		array2d &result
	) const
{
	result.resize( nRows, array1d(nColumns) );

	for (int i=0; i<nRows; i++)
	{
		for (int j=0; j<nColumns; j++)
		{
			result[i][j] = data[i][j] / scaler[i][j];
		}
	}

	return 1;
}

bool toDLib(const ofPixels& inPix, array2d<rgb_pixel>& outPix){
    
    int width = inPix.getWidth();
    int height = inPix.getHeight();
    outPix.set_size( height, width );
    int chans = inPix.getNumChannels();
    const unsigned char* data = inPix.getData();
    for ( unsigned n = 0; n < height;n++ )
    {
        const unsigned char* v =  &data[n * width *  chans];
        for ( unsigned m = 0; m < width;m++ )
        {
            if ( chans==1 )
            {
                unsigned char p = v[m];
                assign_pixel( outPix[n][m], p );
            }
            else{
                rgb_pixel p;
                p.red = v[m*3];
                p.green = v[m*3+1];
                p.blue = v[m*3+2];
                assign_pixel( outPix[n][m], p );
            }
        }
    }
//    if(inPix.getNumChannels() == 3){
//        int h = inPix.getHeight();
//        int w = inPix.getWidth();
//        outPix.clear();
//        outPix.set_size(h,w);
//        for (int i = 0; i < h; i++) {
//            for (int j = 0; j < w; j++) {
//                
//                outPix[i][j].red = inPix.getColor(j, i).r; //inPix[i*w + j];
//                outPix[i][j].green = inPix.getColor(j, i).g; //inPix[i*w + j + 1];
//                outPix[i][j].blue = inPix.getColor(j, i).b; //inPix[i*w + j + 2];
//            }
//        }
//        return true;
//    }else{
//        return  false;
//    }
    return true;
}

static int d8_FlowDir(const array2d<T> &elevations, const int x, const int y){
  T minimum_elevation=elevations(x,y);
  int flowdir=NO_FLOW;

  if (elevations.edge_grid(x,y)){
    if(x==0 && y==0)
      return 2;
    else if(x==0 && y==elevations.height()-1)
      return 8;
    else if(x==elevations.width()-1 && y==0)
      return 4;
    else if(x==elevations.width()-1 && y==elevations.height()-1)
      return 6;
    else if(x==0)
      return 1;
    else if(x==elevations.width()-1)
      return 5;
    else if(y==0)
      return 3;
    else if(y==elevations.height()-1)
      return 7;
  }

  /*NOTE: Since the very edges of the DEM are defined to always flow outwards,
  if they have defined elevations, it is not necessary to check if a neighbour
  is IN_GRID in the following
  NOTE: It is assumed that the no_data datum is an extremely negative
  number, such that all water which makes it to the edge of the DEM's region
  of defined elevations is sucked directly off the grid, rather than piling up
  on the edges.*/
  for(int n=1;n<=8;n++)
    if(
      elevations(x+dx[n],y+dy[n])<minimum_elevation
      || (elevations(x+dx[n],y+dy[n])==minimum_elevation
            && flowdir>0 && flowdir%2==0 && n%2==1)
    ){
      minimum_elevation=elevations(x+dx[n],y+dy[n]);
      flowdir=n;
    }

  return flowdir;
}

int output_ascii_data(
  const std::string filename,
  const array2d<T> &output_grid,
  int precision=8
){
  std::ofstream fout;
  std::string outputsep=" ";
  int output_type=OUTPUT_DEM;
  Timer write_time;
  ProgressBar progress;

  write_time.start();

  diagnostic_arg("Opening ASCII output file \"%s\"...",filename.c_str());
  fout.open(filename.c_str());
  if(!fout.is_open()){
    diagnostic("failed!\n");
    exit(-1);  //TODO: Need to make this safer! Don't just close after all that work!
  }
  diagnostic("succeeded.\n");

  //OmniGlyph output
  if(filename.substr(filename.length()-4)==".omg"){
    outputsep="|";
    output_type=OUTPUT_OMG;
    diagnostic("Writing OmniGlyph file header...");
    fout<<"Contents: Pixel array"<<std::endl;
    fout<<std::endl;
    fout<<"Width:    "<<output_grid.width()<<std::endl;
    fout<<"Height:   "<<output_grid.height()<<std::endl;
    fout<<std::endl;
    fout<<"Spectral bands:   1"<<std::endl;
    fout<<"Bits per band:   32"<<std::endl;
    fout<<"Range of values:   "<<output_grid.min()<<","<<output_grid.max()<<std::endl;
    fout<<"Actual range:   "<<output_grid.no_data<<","<<output_grid.max()<<std::endl;  //TODO: Assumes no_data is a small negative value
    fout<<"Gamma exponent:   0."<<std::endl;
    fout<<"Resolution:   100 pixels per inch"<<std::endl;
    fout<<std::endl;
    fout<<"|"<<std::endl;
  } else {
    diagnostic("Writing ArcGrid ASCII file header...");
    fout<<"ncols\t\t"<<output_grid.width()<<std::endl;
    fout<<"nrows\t\t"<<output_grid.height()<<std::endl;
    fout<<"xllcorner\t"<<std::fixed<<std::setprecision(precision)<<output_grid.xllcorner<<std::endl;
    fout<<"yllcorner\t"<<std::fixed<<std::setprecision(precision)<<output_grid.yllcorner<<std::endl;
    fout<<"cellsize\t"<<std::fixed<<std::setprecision(precision)<<output_grid.cellsize<<std::endl;
    fout<<"NODATA_value\t"<<std::fixed<<std::setprecision(precision);
    if(sizeof(T)==1)  //TODO: Crude way of detecting chars and bools
      fout<<(int)output_grid.no_data<<std::endl;
    else
      fout<<output_grid.no_data<<std::endl;
  }
  diagnostic("succeeded.\n");

  diagnostic("%%Writing ArcGrid ASCII file data...\n");
  progress.start( output_grid.width()*output_grid.height() );
  fout.precision(precision);
  fout.setf(std::ios::fixed);
  for(int y=0;y<output_grid.height();y++){
    progress.update( y*output_grid.width() );
    if(output_type==OUTPUT_OMG)
      fout<<"|";
    for(int x=0;x<output_grid.width();x++)
      if(sizeof(T)==1)  //TODO: Crude way of detecting chars and bools
        fout<<(int)output_grid(x,y)<<outputsep;
      else
        fout<<output_grid(x,y)<<outputsep;
    fout<<std::endl;
  }
  diagnostic_arg(SUCCEEDED_IN,progress.stop());

  fout.close();

  write_time.stop();
  diagnostic_arg("Write time was: %lf\n", write_time.accumulated());

  return 0;
}

namespace amorphous
{

using boost::shared_ptr;


class CStageLightAttributeHolder
{
	enum DefaultBaseEntityHandleIndex
	{
		DIRECTIONAL_LIGHT,
		POINT_LIGHT,
		SPOTLIGHT,
		HS_DIRECTIONAL_LIGHT,
		HS_POINT_LIGHT,
		NUM_DEFAULT_BASE_ENTITY_HANDLES
	};

public:

	boost::weak_ptr<CStage> m_pStage;

	BaseEntityHandle m_aBaseEntityHandle[NUM_DEFAULT_BASE_ENTITY_HANDLES];

	/// target light entity
	EntityHandle<LightEntity> m_TargetLightEntity;
};


CStageLightAttributeHolder gs_StageLightAttribute;


CCopyEntity *CreateEntityFromDesc( CCopyEntityDesc& desc )
{
	CStage *pStage = GetStageForScriptCallback();

	if( pStage )
		return pStage->CreateEntity( desc );
	else
		return NULL;
};



/// row:    holds light groups
/// column: holds object groups (lighting groups)
array2d<char> m_LightToObject;

int init_light_groups()
{
	const int num_light_groups = 16;
	const int num_ligting_groups = 16;
	m_LightToObject.resize( num_light_groups, num_ligting_groups );

	return 0;
}


inline static shared_ptr<LightEntity> GetTargetLightEntity()
{
	return gs_StageLightAttribute.m_TargetLightEntity.Get();
}


inline static void GetInvalidDesc( LightEntityDesc& desc )
{
	// set all colors to invalid values
	desc.aColor[0] = CBE_Light::ms_InvalidColor;
	desc.aColor[1] = CBE_Light::ms_InvalidColor;
	desc.aColor[2] = CBE_Light::ms_InvalidColor;
	desc.LightGroup = CBE_Light::ms_InvalidLightGroup;

	// attenuation (for point lights)
	desc.afAttenuation[0] = desc.afAttenuation[0];
	desc.afAttenuation[1] = desc.afAttenuation[1];
	desc.afAttenuation[2] = desc.afAttenuation[2];
}



using namespace py::light;


static CCopyEntity *GetEntityByName( const char* entity_name )
{
	if( GetStageForScriptCallback() )
        return GetStageForScriptCallback()->GetEntitySet()->GetEntityByName(entity_name);
	else
        return NULL;
}


PyObject* py::light::CreateDirectionalLight( PyObject* self, PyObject* args, PyObject *keywords )
{
	LightEntityDesc desc( Light::DIRECTIONAL );
	char *base_name = "";
	char *light_name = "";
	int shadow_for_light = 1; // true(1) by default
	Vector3 dir = Vector3(0,-1,0); // default direction = vertically down
	SFloatRGBAColor color = SFloatRGBAColor(1,1,1,1);
//.........这里部分代码省略.........

void d8_upslope_cells(
  int x0, int y0, int x1, int y1,
  const array2d<T> &flowdirs,array2d<U> &upslope_cells
){
  diagnostic("Setting up the upslope_cells matrix...");
  upslope_cells.copyprops(flowdirs);
  upslope_cells.init(d8_NO_DATA);
  upslope_cells.no_data=d8_NO_DATA;
  diagnostic("succeeded.\n");
  ProgressBar progress;

  std::queue<grid_cell> expansion;

  if(x0>x1){
    std::swap(x0,x1);
    std::swap(y0,y1);
  }

  //Modified Bresenham Line-Drawing Algorithm
  int deltax=x1-x0;
  int deltay=y1-y0;
  float error=0;
  float deltaerr=(float)deltay/(float)deltax;
  if (deltaerr<0)
    deltaerr=-deltaerr;
  diagnostic_arg("Line slope is %f\n",deltaerr);
  int y=y0;
  for(int x=x0;x<=x1;x++){
    expansion.push(grid_cell(x,y));
    upslope_cells(x,y)=2;
    error+=deltaerr;
    if (error>=0.5) {
      expansion.push(grid_cell(x+1,y));
      upslope_cells(x+1,y)=2;
      y+=sgn(deltay);
      error-=1;
    }
  }

  progress.start(flowdirs.data_cells);
  long int ccount=0;
  while(expansion.size()>0){
    grid_cell c=expansion.front();
    expansion.pop();

    progress.update(ccount++);

    for(int n=1;n<=8;n++)
      if(!flowdirs.in_grid(c.x+dx[n],c.y+dy[n]))
        continue;
      else if(flowdirs(c.x+dx[n],c.y+dy[n])==NO_FLOW)
        continue;
      else if(flowdirs(c.x+dx[n],c.y+dy[n])==flowdirs.no_data)
        continue;
      else if(upslope_cells(c.x+dx[n],c.y+dy[n])==upslope_cells.no_data && n==inverse_flow[flowdirs(c.x+dx[n],c.y+dy[n])]){
        expansion.push(grid_cell(c.x+dx[n],c.y+dy[n]));
        upslope_cells(c.x+dx[n],c.y+dy[n])=1;
      }
  }
  diagnostic_arg(SUCCEEDED_IN,progress.stop());
  diagnostic_arg("Found %ld up-slope cells.\n",ccount);
}

int write_floating_data(
  const std::string basename,
  const array2d<T> &output_grid
){
  Timer write_time;
  ProgressBar progress;
  std::string fn_header(basename), fn_data(basename);

  //TODO: The section below should work, but is something of an abomination
  if(typeid(T)==typeid(float)){
    fn_header+=".hdr";
    fn_data+=".flt";
  } else if (typeid(T)==typeid(double)){
    fn_header+=".hdr";
    fn_data+=".dflt";
  } else {
    std::cerr<<"Cannot write floating type data into this format!"<<std::endl;
    exit(-1);
  }

  write_time.start();


  {
    diagnostic_arg("Opening floating-point header file \"%s\" for writing...",fn_header.c_str());
    std::ofstream fout;
    fout.open(fn_header.c_str());
    if(!fout.is_open()){
      diagnostic("failed!\n");
      exit(-1);  //TODO: Need to make this safer! Don't just close after all that work!
    }
    diagnostic("succeeded.\n");

    diagnostic("Writing floating-point header file...");
    fout<<"ncols\t\t"<<output_grid.width()<<std::endl;
    fout<<"nrows\t\t"<<output_grid.height()<<std::endl;
    fout<<"xllcorner\t"<<std::fixed<<std::setprecision(10)<<output_grid.xllcorner<<std::endl;
    fout<<"yllcorner\t"<<std::fixed<<std::setprecision(10)<<output_grid.yllcorner<<std::endl;
    fout<<"cellsize\t"<<std::fixed<<std::setprecision(10)<<output_grid.cellsize<<std::endl;
    fout<<"NODATA_value\t"<<std::fixed<<std::setprecision(10)<<output_grid.no_data<<std::endl;
    fout<<"BYTEORDER\tLSBFIRST"<<std::endl; //TODO
    fout.close();
    diagnostic("succeeded.\n");
  }


  diagnostic_arg("Opening floating-point data file \"%s\" for writing...",fn_data.c_str());

  {
    std::ofstream fout(fn_data.c_str(), std::ios::binary | std::ios::out);
    if(!fout.is_open()){
      diagnostic("failed!\n");
      exit(-1);  //TODO: Need to make this safer! Don't just close after all that work!
    }
    diagnostic("succeeded.\n");

    diagnostic("%%Writing floating-point data file...\n");
    progress.start( output_grid.width()*output_grid.height() );
    for(int y=0;y<output_grid.height();++y){
      progress.update( y*output_grid.width() );
      for(int x=0;x<output_grid.width();++x)
        fout.write(reinterpret_cast<const char*>(&output_grid(x,y)), std::streamsize(sizeof(T)));
    }
    fout.close();
    write_time.stop();
    diagnostic_arg(SUCCEEDED_IN,progress.stop());
  }

  diagnostic_arg("Write time was: %lf\n", write_time.accumulated());

  return 0;
}