function [HME]=hmegen(params,k)
	#Params=[C0;-T2;-T1;Omega;Omegax;Omegay]
	spins=columns(params);
	size(params)
	C0=params(1,:);
	T1=params(3,:);
	T2=params(2,:);
	Omega=params(4,:);
	Omegax=params(5,:);
	Omegay=params(6,:);
	
	for n=1:spins
		M{n}=diag([-1/T2(n) -1/T2(n) -1/T2(n)])+diag([Omega(n) Omegax(n)],1)-diag([Omega(n) Omegax(n)],-1)+diag([-Omegay(n)],2)+diag(Omegay(n),-2);
	endfor
	

	HME=zeros(3*n,3*n);
	size(HME);
	for n=1:spins

			HME((n-1)*3+1:(n-1)*3+3,(n-1)*3+1:(n-1)*3+3)=M{n};
	endfor
	HME=HME+k;
	HME=[zeros(1,columns(HME)); HME];

	phi=[0];
	for n=1:spins
		phi=[phi 0 0 C0(n)/(-T1(n))];
	endfor
	spins;
	phi';
	size(HME);
	HME=[-2.*phi' HME]
endfunction
