function fterm=add_terms(Sb,So,ii,yy,rindice,op,model,H,dir)

Ni=length(rindice);
Nc=set_numberofcombinations(rindice,model);
stindex=set_startindex(yy,dir);

tmpadjS=Sb(rindice,:);
tmpadjS=repmat(tmpadjS,[2,1]);
tmpadjS(Ni+1:end,:)=-tmpadjS(Ni+1:end,:);

if(strcmp(dir,'horizontal')==1)
    index=(yy-1)*op.xdim+ii+stindex-1;
    Snei=extract_neigboring(So,index,H);
    
    if(yy==1)
        if(ii==1)    
            if(strcmp(model,'TICA')==1)
                tmpadjS=Sb(rindice,:);
            end
            
            tmpS=repmat(Sb(1,:),[Nc,1]);
            fterm=compute_comb(tmpS,tmpadjS,model);
        else
            fterm=0;
        end
    elseif ((1<yy)&&(yy<op.ydim))   
        if(op.toponum==8)
            % upper left component
            term1=compute_comb(repmat(Snei(1,:),[Nc,1]),tmpadjS,model);        
            % upper component
            term2=compute_comb(repmat(Snei(2,:),[Nc,1]),tmpadjS,model);        
            % upper right component
            term3=compute_comb(repmat(Snei(3,:),[Nc,1]),tmpadjS,model);
        
            fterm=term1+term2+term3;

            if(ii==1)
                % left component
                term4=compute_comb(repmat(Snei(4,:),[Nc,1]),tmpadjS,model);        
                % below left component
                term5=compute_comb(repmat(Snei(6,:),[Nc,1]),tmpadjS,model);
            
                fterm=fterm+term4+term5;
            end
        else
            % upper component
            term1=compute_comb(repmat(Snei(1,:),[Nc,1]),tmpadjS,model);        
            
            fterm=term1;
            
            if(ii==1)
                % left component
                term2=compute_comb(repmat(Snei(2,:),[Nc,1]),tmpadjS,model);        
            
                fterm=fterm+term2;
            end
        end
    else       
        if(op.toponum==8)
            % upper left component
            term1=compute_comb(repmat(Snei(1,:),[Nc,1]),tmpadjS,model);        
            % upper component
            term2=compute_comb(repmat(Snei(2,:),[Nc,1]),tmpadjS,model);        
            % upper right component
            term3=compute_comb(repmat(Snei(3,:),[Nc,1]),tmpadjS,model);
            % below left component
            term4=compute_comb(repmat(Snei(6,:),[Nc,1]),tmpadjS,model);
            % below component
            term5=compute_comb(repmat(Snei(7,:),[Nc,1]),tmpadjS,model);        
            % below right component
            term6=compute_comb(repmat(Snei(8,:),[Nc,1]),tmpadjS,model);

            fterm=term1+term2+term3+term4+term5+term6;
        
            if(ii==1)
                % left component
                term7=compute_comb(repmat(Snei(4,:),[Nc,1]),tmpadjS,model);
                    
                fterm=fterm+term7;
            end
        else
            % upper component
            term1=compute_comb(repmat(Snei(1,:),[Nc,1]),tmpadjS,model);        
            % below component
            term2=compute_comb(repmat(Snei(4,:),[Nc,1]),tmpadjS,model);        

            fterm=term1+term2;
            
            if(ii==1)
                % left component
                term3=compute_comb(repmat(Snei(2,:),[Nc,1]),tmpadjS,model);
                
                fterm=fterm+term3;
            end
        end
    end
elseif(strcmp(dir,'vertical')==1)
    index=(ii+stindex-2)*op.xdim+yy;
    Snei=extract_neigboring(So,index,H);
    
    if(yy==1)
        if(ii==1)
            if(op.toponum==8)
                % upper left component
                term1=compute_comb(repmat(Snei(1,:),[Nc,1]),tmpadjS,model);
                % upper component
                term2=compute_comb(repmat(Snei(2,:),[Nc,1]),tmpadjS,model);
                % upper right component
                term3=compute_comb(repmat(Snei(3,:),[Nc,1]),tmpadjS,model);
                
                fterm=term1+term2+term3;
            else
                % upper component
                term1=compute_comb(repmat(Snei(1,:),[Nc,1]),tmpadjS,model);
                
                fterm=term1;
            end                
        else
            fterm=0;            
        end
    else
        if(op.toponum==8)
            % upper left component
            term1=compute_comb(repmat(Snei(1,:),[Nc,1]),tmpadjS,model);
            % left component
            term2=compute_comb(repmat(Snei(4,:),[Nc,1]),tmpadjS,model);
            % below left component
            term3=compute_comb(repmat(Snei(6,:),[Nc,1]),tmpadjS,model);
            
            fterm=term1+term2+term3;

            if(ii==1)
                % upper component
                term4=compute_comb(repmat(Snei(2,:),[Nc,1]),tmpadjS,model);                
                % upper right component
                term5=compute_comb(repmat(Snei(3,:),[Nc,1]),tmpadjS,model);
                
                fterm=fterm+term4+term5;
            end
        else
            % left component
            term1=compute_comb(repmat(Snei(2,:),[Nc,1]),tmpadjS,model);

            fterm=term1;
            
            if(ii==1)
                % upper component
                term2=compute_comb(repmat(Snei(1,:),[Nc,1]),tmpadjS,model);
                
                fterm=fterm+term2;
            end
        end
    end
else
    error('Specify the direction you want to optimize.\n');
end