%MESH2DVERTICAL class definition
%
%   Usage:
%      mesh2dvertical=mesh2dvertical();

classdef mesh2dvertical
	properties (SetAccess=public) 
		x                           = NaN;
		y                           = NaN;
		elements                    = NaN
		numberofelements            = 0;
		numberofvertices            = 0;
		numberofedges               = 0;

		lat                         = NaN
		long                        = NaN
		hemisphere                  = NaN

		vertexonboundary            = NaN
		vertexonbase                 = NaN
		vertexonsurface             = NaN

		edges                       = NaN
		segments                    = NaN
		segmentmarkers              = NaN
		vertexconnectivity          = NaN
		elementconnectivity         = NaN
		average_vertex_connectivity = 0;
	end
	methods
                function createxml(obj,fid) % {{{
            fprintf(fid, '<!-- 2d Mesh (Vertical) -->\n');
            
            %elements and vertices
            fprintf(fid,'%s\n%s\n%s\n','       <frame key="1" label="Elements and vertices">','       <section name="mesh" />');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n','          <parameter key ="numberofelements" type="',class(obj.numberofelements),'" default="',convert2str(obj.numberofelements),'">','              <help> number of elements </help>','          </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n','          <parameter key ="numberofvertices" type="',class(obj.numberofvertices),'" default="',convert2str(obj.numberofvertices),'">','              <help> number of vertices </help>','          </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n','          <parameter key ="elements" type="',class(obj.elements),'" default="',convert2str(obj.elements),'">','              <help> vertex indices of the mesh elements </help>','          </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="x" type="',class(obj.x),'" default="',convert2str(obj.x),'">','     <section name="mesh" />','     <help> vertices x coordinate [m] </help>','  </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="y" type="',class(obj.y),'" default="',convert2str(obj.y),'">','     <section name="mesh" />','     <help> vertices y coordinate [m] </help>','  </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="edges" type="',class(obj.edges),'" default="',convert2str(obj.edges),'">','     <section name="mesh" />','     <help> edges of the 2d mesh (vertex1 vertex2 element1 element2) </help>','  </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="numberofedges" type="',class(obj.numberofedges),'" default="',convert2str(obj.numberofedges),'">','     <section name="mesh" />','     <help> number of edges of the 2d mesh </help>','  </parameter>');
            fprintf(fid,'%s\n%s\n','       </frame>');
     
            % properties
            fprintf(fid,'%s\n%s\n%s\n','<frame key="2" label="Properties">','<section name="mesh" />');             
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="vertexonboundary" type="',class(obj.vertexonboundary),'" default="',convert2str(obj.vertexonboundary),'">','     <section name="mesh" />','     <help> vertices on the boundary of the domain flag list </help>','  </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="vertexonbase" type="',class(obj.vertexonbase),'" default="',convert2str(obj.vertexonbase),'">','     <section name="mesh" />','     <help> vertices on the bed of the domain flag list </help>','  </parameter>');            
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="segments" type="',class(obj.segments),'" default="',convert2str(obj.segments),'">','     <section name="mesh" />','     <help> edges on domain boundary (vertex1 vertex2 element) </help>','  </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="segmentmarkers" type="',class(obj.segmentmarkers),'" default="',convert2str(obj.segmentmarkers),'">','     <section name="mesh" />','     <help> number associated to each segment </help>','  </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="vertexconnectivity" type="',class(obj.vertexconnectivity),'" default="',convert2str(obj.vertexconnectivity),'">','     <section name="mesh" />','     <help> list of vertices connected to vertex_i </help>','  </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="elementconnectivity" type="',class(obj.elementconnectivity),'" default="',convert2str(obj.elementconnectivity),'">','     <section name="mesh" />','     <help> list of vertices connected to element_i </help>','  </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="average_vertex_connectivity" type="',class(obj.average_vertex_connectivity),'" default="',convert2str(obj.average_vertex_connectivity),'">','     <section name="mesh" />','     <help> average number of vertices connected to one vertex </help>','  </parameter>');
            fprintf(fid,'%s\n%s\n','</frame>');

            %projection
            fprintf(fid,'%s\n%s\n%s\n','<frame key="3" label="Projection">','<section name="mesh" />'); 
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="lat" type="',class(obj.lat),'" default="',convert2str(obj.lat),'">','     <section name="mesh" />','     <help> vertices latitude [degrees] </help>','  </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="long" type="',class(obj.long),'" default="',convert2str(obj.long),'">','     <section name="mesh" />','     <help> verticies longitude [degrees] </help>','  </parameter>');
            fprintf(fid,'%s%s%s%s%s\n%s\n%s\n','  <parameter key ="hemisphere" type="',class(obj.hemisphere),'" default="',convert2str(obj.hemisphere),'">','     <section name="mesh" />','     <help> Indicate hemisphere ''n'' or ''s'' </help>','  </parameter>');
            fprintf(fid,'%s\n%s\n','</frame>');
        
        end % }}}
		function obj = mesh2dvertical(varargin) % {{{
			switch nargin
				case 0
					obj=setdefaultparameters(obj);
				case 1
					inputstruct=varargin{1};
					list1 = properties('mesh2dvertical');
					list2 = fieldnames(inputstruct);
					for i=1:length(list1)
						fieldname = list1{i};
						if ismember(fieldname,list2),
							obj.(fieldname) = inputstruct.(fieldname);
						end
					end
				otherwise
					error('constructor not supported');
			end
		end % }}}
		function obj = setdefaultparameters(obj) % {{{

			%the connectivity is the averaged number of nodes linked to a
			%given node through an edge. This connectivity is used to initially
			%allocate memory to the stiffness matrix. A value of 16 seems to
			%give a good memory/time ration. This value can be checked in
			%trunk/test/Miscellaneous/runme.m
			obj.average_vertex_connectivity=25;
		end % }}}
		function md = checkconsistency(obj,md,solution,analyses) % {{{

			md = checkfield(md,'fieldname','mesh.x','NaN',1,'size',[md.mesh.numberofvertices 1]);
			md = checkfield(md,'fieldname','mesh.y','NaN',1,'size',[md.mesh.numberofvertices 1]);
			md = checkfield(md,'fieldname','mesh.elements','NaN',1,'>',0,'values',1:md.mesh.numberofvertices);
			md = checkfield(md,'fieldname','mesh.elements','size',[md.mesh.numberofelements 3]);
			if any(~ismember(1:md.mesh.numberofvertices,sort(unique(md.mesh.elements(:)))));
				md = checkmessage(md,'orphan nodes have been found. Check the mesh outline');
			end
			md = checkfield(md,'fieldname','mesh.numberofelements','>',0);
			md = checkfield(md,'fieldname','mesh.numberofvertices','>',0);
			md = checkfield(md,'fieldname','mesh.vertexonbase','size',[md.mesh.numberofvertices 1],'values',[0 1]);
			md = checkfield(md,'fieldname','mesh.vertexonsurface','size',[md.mesh.numberofvertices 1],'values',[0 1]);
			md = checkfield(md,'fieldname','mesh.average_vertex_connectivity','>=',9,'message','''mesh.average_vertex_connectivity'' should be at least 9 in 2d');

			switch(solution),
				case ThermalSolutionEnum(),
					md = checkmessage(md,'thermal not supported for 2d mesh');
			end
		end % }}}
		function disp(obj) % {{{
			disp(sprintf('   2D tria Mesh (vertical):')); 

			disp(sprintf('\n      Elements and vertices:'));
			fielddisplay(obj,'numberofelements','number of elements');
			fielddisplay(obj,'numberofvertices','number of vertices');
			fielddisplay(obj,'elements','vertex indices of the mesh elements');
			fielddisplay(obj,'x','vertices x coordinate [m]');
			fielddisplay(obj,'y','vertices y coordinate [m]');
			fielddisplay(obj,'edges','edges of the 2d mesh (vertex1 vertex2 element1 element2)');
			fielddisplay(obj,'numberofedges','number of edges of the 2d mesh');

			disp(sprintf('\n      Properties:'));
			fielddisplay(obj,'vertexonboundary','vertices on the boundary of the domain flag list');
			fielddisplay(obj,'vertexonbase','vertices on the bed of the domain flag list');
			fielddisplay(obj,'vertexonsurface','vertices on the surface of the domain flag list');
			fielddisplay(obj,'segments','edges on domain boundary (vertex1 vertex2 element)');
			fielddisplay(obj,'segmentmarkers','number associated to each segment');
			fielddisplay(obj,'vertexconnectivity','list of vertices connected to vertex_i');
			fielddisplay(obj,'elementconnectivity','list of vertices connected to element_i');
			fielddisplay(obj,'average_vertex_connectivity','average number of vertices connected to one vertex');

			disp(sprintf('\n      Projection:'));
			fielddisplay(obj,'lat','vertices latitude [degrees]');
			fielddisplay(obj,'long','vertices longitude [degrees]');
			fielddisplay(obj,'hemisphere','Indicate hemisphere ''n'' or ''s'' ');
		end % }}}
		function marshall(obj,md,fid) % {{{
			WriteData(fid,'enum',DomainTypeEnum(),'data',StringToEnum(['Domain' domaintype(obj)]),'format','Integer');
			WriteData(fid,'enum',DomainDimensionEnum(),'data',dimension(obj),'format','Integer');
			WriteData(fid,'enum',MeshElementtypeEnum(),'data',StringToEnum(elementtype(obj)),'format','Integer');
			WriteData(fid,'object',obj,'class','mesh','fieldname','x','format','DoubleMat','mattype',1);
			WriteData(fid,'object',obj,'class','mesh','fieldname','y','format','DoubleMat','mattype',1);
			WriteData(fid,'enum',MeshZEnum(),'data',zeros(obj.numberofvertices,1),'format','DoubleMat','mattype',1);
			WriteData(fid,'object',obj,'class','mesh','fieldname','elements','format','DoubleMat','mattype',2);
			WriteData(fid,'object',obj,'class','mesh','fieldname','numberofelements','format','Integer');
			WriteData(fid,'object',obj,'class','mesh','fieldname','numberofvertices','format','Integer');
			WriteData(fid,'object',obj,'class','mesh','fieldname','vertexonbase','format','BooleanMat','mattype',1);
			WriteData(fid,'object',obj,'class','mesh','fieldname','vertexonsurface','format','BooleanMat','mattype',1);
			WriteData(fid,'object',obj,'class','mesh','fieldname','average_vertex_connectivity','format','Integer');
		end % }}}
		function t = domaintype(obj) % {{{
			t = '2Dvertical';
		end % }}}
		function d = dimension(obj) % {{{
			d = 2;
		end % }}}
		function s = elementtype(obj) % {{{
			s = 'Tria';
		end % }}}
		function flags = vertexflags(self,value) % {{{
			flags = zeros(self.numberofvertices,1);
			pos   = self.segments(find(self.segmentmarkers==value),1:2);
			flags(pos) = 1;
		end % }}}
		function [data datatype] = processdata(self,md,data,options) % {{{

			%transpose data if necessary
			if (size(data,2) > size(data,1)),
				data=data';
			end
			datasize=size(data);

			%convert to double if necessary
			if ~isnumeric(data);
				disp('processdata info message: data is not numeric (logical?). Converted to double');
				data=double(data);
			end

			%check length
			if datasize(1)~=md.mesh.numberofvertices & datasize(1)~=md.mesh.numberofelements
				error('plotmodel error message: data not supported yet');
			end

			%quiver?
			if datasize(2)>1,
				datatype=3;
			end

			%smoothing?
			if exist(options,'smooth')
				data=averaging(md,data,getfieldvalue(options,'smooth'));
				datasize(1)=md.mesh.numberofvertices;
				%---> go to node data
			end

			%element data
			if (datasize(1)==md.mesh.numberofelements & datasize(2)==1),
				datatype=1;

				%Mask?
				if exist(options,'mask'),
					flags=getfieldvalue(options,'mask');
					pos=find(~flags);
					if length(flags)==md.mesh.numberofvertices,
						[pos2 dummy]=find(ismember(md.mesh.elements,pos));
						data(pos2,:)=NaN;
					elseif length(flags)==md.mesh.numberofelements
						data(pos,:)=NaN;
					else
						disp('plotmodel warning: mask length not supported yet (supported length are md.mesh.numberofvertices and md.mesh.numberofelements');
					end
				end

				%log?
				if exist(options,'log'),
					bounds=getfieldvalue(options,'caxis',[min(data(:)) max(data(:))]);
					data(find(data<bounds(1)))=bounds(1);
					if any(data<=0),
						error('Log option cannot be applied on negative values. Use caxis option (Rignot''s settings: [1.5 max(data)])');
					end
					pos=find(~isnan(data));
					data(pos)=log(data(pos))/log(getfieldvalue(options,'log'));
				end
			end

			%node data
			if (datasize(1)==md.mesh.numberofvertices & datasize(2)==1),
				datatype=2;

				%Mask?
				if exist(options,'mask'),
					flags=getfieldvalue(options,'mask');
					pos=find(~flags);
					if length(flags)==md.mesh.numberofvertices,
						data(pos,:)=NaN;
					elseif length(flags)==md.mesh.numberofelements
						data(md.mesh.elements(pos,:),:)=NaN;
					else
						disp('plotmodel warning: mask length not supported yet (supported length are md.mesh.numberofvertices and md.mesh.numberofelements');
					end
				end

				%log?
				if exist(options,'log'),
					%if any(data<=0),
					%	error('Log option cannot be applied on negative values. Use caxis option (Rignot''s settings: [1.5 max(data)])');
					%end
					data=log(data)/log(getfieldvalue(options,'log'));
				end
			end
		end % }}}
		function [x y z elements is2d isplanet] = processmesh(self,options) % {{{

			isplanet = 0;
			is2d     = 1;

			elements = self.elements;
			x        = self.x;
			y        = self.y;
			z        = zeros(self.numberofvertices,1);
		end % }}}
	end
end
