import numpy as  np

def project2d(md3d,value,layer):
	'''
	returns the value of a field for a given layer of the mesh
	

   returns the value of a vector for a given layer from extruded mesh onto the 2d mesh 
   used to do the extrusion. This function is used to compare values between different
   layers of a 3d mesh.

   Usage:
      projection_value=project2d(md3d,value,layer)

   Example:
      vel2=project2d(md3d,md3d.initialization.vel,2);
      returns the velocity of the second layer (1 is the base)
	'''

	if md3d.mesh.domaintype().lower() != '3d':
		raise StandardError("model passed to project2d function should be 3D")

	if layer<1 or layer>md3d.mesh.numberoflayers:
		raise ValueError("layer must be between 0 and %i" % md3d.mesh.numberoflayers)
	
	# coerce to array in case float is passed
	if type(value)!=np.ndarray:
		print 'coercing array'
		value=np.array(value)

	vec2d=False
	if value.ndim==2 and value.shape[1]==1: 
		value=value.reshape(-1,)
		vec2d=True

	if value.size==1:
		projection_value=value[(layer-1)*md3d.mesh.numberofelements2d:layer*md3d.mesh.numberofelements2d]
	elif value.shape[0]==md3d.mesh.numberofvertices:
		#print 'indices: ', (layer-1)*md3d.mesh.numberofvertices2d, layer*md3d.mesh.numberofvertices2d
		projection_value=value[(layer-1)*md3d.mesh.numberofvertices2d:layer*md3d.mesh.numberofvertices2d]
	elif value.shape[0]==md3d.mesh.numberofvertices+1:
		projection_value=[value[(layer-1)*md3d.mesh.numberofvertices2d:layer*md3d.mesh.numberofvertices2d], value[-1]]
	else:
		projection_value=value[(layer-1)*md3d.mesh.numberofelements2d:layer*md3d.mesh.numberofelements2d]

	if vec2d:
		projection_value=projection_value.reshape(-1,)

	return projection_value
