function MineralStruct = CreateMineralCell(Minerals,NumberMinerals,thelabels,controlstruct,which0,MineralStructIn,InStruct)
%function MineralCell = CreateMineralCell(Minerals,NumberMinerals,thelabels,controlstruct,which0)
%This function creates a cell array where the columns are the minerals that will be considered
%while the rows are the data that is stored.  The neccessary vectors and matrices are shortened
%using the function shorton depending on the actual elements present in the input database
%
%INPUT
%Minerals = A cell array of the names of the minerals that are being considered for identification
%NumberMinerals = The number of minerals that are considered
%thelabels = The list of elements that are included in database
%controlstruct = Contains options info for this run (created by APHID.m)
%which0 = Information about which oxygen value to use (hydrous or anhydrous)
%MineralStructIn = The current state of MineralStuct, this function modifies MineralStructIn
%and then outputs it as MineralStruct
%InStuct = the raw data read from the background file, created by loadbackground.m  This function
%uses InStruct as background information to modify MineralStruct.
%
%OUTPUT
%MineralStruct = a structure array that contains all of the information required
%to convert a weights to cations and finally to additive/exchange and then back
%
%Written by K. Steffen (ksteffen@unm.edu) at the University of New Mexico,
%Department of Earch and Planetary Science
%Copyright 2003 Kurt Steffen
%
%Last Modified June 17, 2003

    
thelabelsarchive = thelabels;
for k = 1:NumberMinerals
	thelabels = thelabelsarchive;  %initialization phase
	irontemp = MineralStructIn(k).IronOxidation;
	irontemp = num2str(irontemp);
	A = MineralStructIn(k).A;
	C = MineralStructIn(k).C;
	D = MineralStructIn(k).D;
	AdditiveExchangeVector = MineralStructIn(k).AdditiveExchangeVector;
	Limits = [MineralStructIn(k).upperlimit ; MineralStructIn(k).lowerlimit];
	order = MineralStructIn(k).order;
	EndMembers = MineralStructIn(k).EndMemberVector;
	Endbase = MineralStructIn(k).EndBase;
 	Conversion  = MineralStructIn(k).EndMemberConversion;
      
	EndMemberCell(2,k) = {EndMembers};
 	EndMemberCell(3,k) = {Conversion};
	[Anew Dnew neworder newLimits newEndMembers newConversion newEndbase MineralStructIn(k).EndMemberCations,AdditiveExchangeVector] = ...
	shorten(A,C,D,order,thelabels,Limits,EndMembers,Conversion,Endbase,MineralStructIn(k).EndMemberCations,AdditiveExchangeVector);
    
	if isempty(neworder) %The components present in the endmember must all also be in the system
		error(['Given the components, the phase ' [Minerals{k}] ' can not be found in system']);
	else
		newordertemp = neworder(1,1:end-MineralStructIn(k).NumberAdditiveComponents);
		mapping = GetColumnsGeneral(thelabels,newordertemp);
		mapping = [mapping{2,:}];
    
		temp = zeros(1,MineralStructIn(k).NumberEndMembers);
		for m =1:length(newEndMembers)
			for n = 1:MineralStructIn(k).NumberEndMembers
				if strcmp(newEndMembers{m},MineralStructIn(k).EndMemberVector{n})
					temp(n) = 1;
				end
			end       
		end
   
		MineralStructIn(k).NumberEndMembers = length(newEndMembers);	
		MineralStructIn(k).EndMemberVector = newEndMembers;
		temporder = neworder(1,1:end-MineralStructIn(k).NumberAdditiveComponents);
		templabels = neworder(1,1:end-MineralStructIn(k).NumberAdditiveComponents);
		MineralStructIn(k).NumberCationComponents = length(temporder);


		%Code to check if endmembers are stated correctly
		whichbase = find(newEndbase == 1);
		ThisNumberAdditive = MineralStructIn(k).NumberAdditiveComponents;
		ThisA = Anew(1:end-ThisNumberAdditive,ThisNumberAdditive+1:end);
		newConversiontemp = newConversion(ThisNumberAdditive+1:end,:);
		newConversiontemp = [newConversiontemp(1:whichbase-1,:) ;newEndbase ; newConversiontemp(whichbase:end,:)];
		ConversionInv = inv(newConversiontemp(1:size(newConversiontemp,1),1:size(newConversiontemp,1)));
		whichbase = find(newEndbase == 1);
		base = MineralStructIn(k).EndMemberCations(whichbase,:);
		for m = 1:size(ConversionInv,1)
			thisConversion = ConversionInv(m,:)';
			thisConversion(whichbase) = [];

			if isempty(thisConversion)
				if max(abs(base  - MineralStructIn(k).EndMemberCations(m,:)))  > .01

					disp([{'From Exchange'} {'From Endmember'}  {' '} ;double2cell([(base + (ThisA*thisConversion)')' ...
					MineralStructIn(k).EndMemberCations(m,:)']) newordertemp'])
			
					error(sprintf('In Mineral %s there is a possible error in the definition of endmember %s. See difference above'...
					,MineralStructIn(k).MineralName,newEndMembers{m}))
				end
			else
				if max(abs(base + (ThisA*thisConversion)' - MineralStructIn(k).EndMemberCations(m,:)))  > .01

					disp([{'From Exchange'} {'From Endmember'}  {' '} ;double2cell([(base + (ThisA*thisConversion)')' ...
					MineralStructIn(k).EndMemberCations(m,:)']) newordertemp'])
			
					error(sprintf('In Mineral %s there is a possible error in the definition of endmember %s. See difference above'...
					,MineralStructIn(k).MineralName,newEndMembers{m}))
				end

			end
		end
		%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

		ironcolumn = GetColumnsGeneral3(templabels,{'Fe'});
		ironcolumn = ironcolumn{2,1};
		gflabels = InStruct.GfElement;
		gf = [InStruct.GramFormulaList';InStruct.CationsList';InStruct.ChargeList']; 
		Oxygens = MineralStructIn(k).NumberOxygenHyd;
 		FeOx = MineralStructIn(k).IronOxidation;
		thenumbers = zeros(length(MineralStructIn(k).NumberEndMembers),MineralStructIn(k).NumberCationComponents);
		Molecular = zeros(length(MineralStructIn(k).NumberEndMembers),1);
    
		%determining percentage oxide weights, molecular weights, for each endmember
		for m = 1:MineralStructIn(k).NumberEndMembers
			templabelstemp = temporder;
			if ironcolumn ~= -1
				templabelstemp(ironcolumn) = {['Fe' num2str([FeOx])]};
			end
			out = prepweights2cations(templabelstemp,gflabels,gf);
			out.Oxygens = Oxygens;
        
			Molecular(m) = sum((MineralStructIn(k).EndMemberCations(m,:)./out.cation).*out.gramformula); % the molecular weight of each endmember
			MineralStructIn(k).EndMemberCations(m,:);
			thenumbers(m,:) = cations2weights(out,MineralStructIn(k).EndMemberCations(m,:),100);
		end
		[rows cols] = size(thenumbers);

		if k == 1
			writecsv([pwd filesep 'Output' filesep 'Endmembers.csv'],[neworder(1,1:cols) {'GramFormula'} {[]}; double2cell(thenumbers) double2cell(Molecular')  newEndMembers'],'w');
		else
			writecsv([pwd filesep 'Output' filesep 'Endmembers.csv'],[neworder(1,1:cols) {'GramFormula'} {[]}; double2cell(thenumbers) double2cell(Molecular') newEndMembers'],'a');
		end

  
		%This code tests to make sure the endmembers are stoichiometric
		[datarow datacol] = size(thenumbers);
		datatest = zeros(datarow,datacol);
    
		%takes the weights that are given and converts back to cations 
		for m = 1:MineralStructIn(k).NumberEndMembers
			templabelstemp = temporder;
			if ironcolumn ~= -1
				templabelstemp(ironcolumn) = {['Fe' num2str([FeOx])]};
			end
			out = prepweights2cations(templabelstemp,gflabels,gf);
			out.Oxygens = Oxygens;
			datatest(m,:) = weights2cations(out,thenumbers(m,:));    
		end
		test = max(max(abs(MineralStructIn(k).EndMemberCations-datatest)));
    
		%if there is a difference in the cation values then print an error message
		if test > .01
 			thetable = double2cell(MineralStructIn(k).EndMemberCations);
			thetable = [MineralStructIn(k).EndMemberVector' thetable];

			disp('Original Values');
  			disp([  {0} temporder; thetable])
        
			thetable = double2cell(MineralStructIn(k).EndMemberCations-datatest);
			thetable = [MineralStructIn(k).EndMemberVector' thetable];

			disp('Error');
			disp([  {0} temporder; thetable])
			error('An endMember in the Endmember table is not stoichiometric: Above is the error table');
		end
    
		%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
		if which0 == 1  %selecting which Oxygen value to use
			O = MineralStructIn(k).NumberOxygen;
		else
			O = MineralStructIn(k).NumberOxygenHyd;
		end

		
		%output to structure array
		MineralStruct(k) = struct(... 
		'density',MineralStructIn(k).density,...
		'MineralName',{Minerals(k)},...
		'A',{Anew},...
		'C',{C},...
		'D',{Dnew},...
		'order',{neworder},...
		'O',{O},...
		'AdditiveExchangeVector',{AdditiveExchangeVector},...
		'NumberAdditiveComponents',{MineralStructIn(k).NumberAdditiveComponents},...
		'NumberExchangeComponents',{MineralStructIn(k).NumberExchangeComponents},...
		'NumberCationComponents',{MineralStructIn(k).NumberCationComponents},...
		'Molecular',{Molecular},...
		'irontemp',{irontemp},...
		'Limits',{newLimits},...
		'EndMembers',{newEndMembers},...
		'Conversion',{newConversion},...
		'EndMemberWeights',{thenumbers},...
		'EndMemberMolecular',{Molecular},...
		'mapping',{mapping},...
		'Endbase',{newEndbase});

	end
end

if NumberMinerals == 0  %If somehow there are no minerals then initialize a Null structure
	MineralStruct = struct([]);
end