checkNetState(netState , phys, trafficMatrix , flag_checkLoops, flag_removeLoops , flag_checkLpRoutingMatrix , flag_lossesAllowed)

% checkNetState
%
% Usage: [flag_modifiedNetState newNetState] = checkNetState(netState , phys, trafficMatrix , ...
%        flag_checkLoops , flag_removeLoops , flag_checkLpRoutingMatrix , flag_lossesAllowed)
%
%The function checks if the 'netState' structure is right and coherent with
% the 'phys' structure and with the 'trafficMatrix'. If there is a error on
% checking, a error is returned, except if the the error is a detected loop
% and the input 'flag_removeLoops'=1, then, the loop is removed.
%
%INPUT arguments:
    %1)netState: original 'netState' structure
    %2)phys: 'phys' structure 
    %3)trafficMatrix: offered traffic matrix. 
    %NOTE: If this field is an empty matrix, the offered traffic is
    %computed from the 'netState.flowTable'. We don't need to use
    %'trafficmatrix' in dynamic scenarios.
    %4)flag_checkLoops: '1' if the function is allowed to check if there
    % are loops in both lightpaths and flows. '0' if the function do not
    % check if there exist loops.
    %5)flag_removeLoops: '1' if the function is allowed to remove loops in
    % both lightpaths and flows. '0', otherwise.
    %6)flag_checkLpRoutingMatrix: '1' if the function is allowed to check 
    % the lightpath routing matrix. '0', otherwise.
    %7)flag_lossesAllowed: '1' if traffic losses are allowed in the flow
    % routing. '0', otherwise.
    
%OUTPUT arguments:
    %1)flag_modifiedNetState: '1', if the 'netState' structure was modified
    %  on removing loops; '0', if 'newNetState' = 'netState'. 
    %2)newNetState: outgoing 'netState' structure. It can be the same as
    %  the input 'netState', or the sctruc without loops if the input option
    %  'flag_checkLpRoutingMatrix' is activated and loops were detected.
    
function [flag_modifiedNetState newNetState]= checkNetState(netState , phys, trafficMatrix , flag_checkLoops, flag_removeLoops , flag_checkLpRoutingMatrix , flag_lossesAllowed)

flag_modifiedNetState = 0; 
newNetState = netState;

%%%% 1. Checking sizes and dimensions in netState %%%%
nrOfNodes = phys.N;
nrOfFibres = phys.M;
nrOfLightpaths= size(netState.lightpathTable,1);
nrOfFlows= size(netState.flowTable,1);

if ((nrOfNodes < 1) ||(nrOfFibres < 1) || (nrOfLightpaths < 1) || (nrOfFlows < 1)), error ('Some field is empty'); end
if not(isequal(size(netState.lightpathTable),[nrOfLightpaths 3])), error('Wrong dimensions size of "lightpathTable"'); end
if not(isequal(size(netState.flowTable),[nrOfFlows 7])), error('Wrong dimensions size of "flowTable"'); end 

if (flag_checkLpRoutingMatrix == 1) && not(isequal(size(netState.lightpathRoutingMatrix),[nrOfLightpaths nrOfFibres])),
    error ('Wrong Size in "lightpathRoutingMatrix"');
end

if not(isequal(size(netState.flowRoutingMatrix),[nrOfFlows nrOfLightpaths])),
    error ('Wrong Size in "flowRoutingMatrix"');
end

if (numel(netState.numberOfOccupiedTxs)~=nrOfNodes || numel(netState.numberOfOccupiedRxs)~=nrOfNodes || numel(netState.numberOfOccupiedTWCs)~=nrOfNodes),
    error ('Wrong Size in "numberOfOccupiedTxs" or in "numberOfOccupiedRxs" or in "numberOfOccupiedTWCs"');
end

%%%% 2. Checking the wavelengths used in the fibre links %%%%
if (flag_checkLpRoutingMatrix ~= 0) 
    for m = 1:nrOfFibres,
        % Checking the Wavelength Validity in the current fibre    
        if sum(ismember(netState.lightpathRoutingMatrix(:,m), 0:phys.numberWavelengthPerFiber(m))) ~= nrOfLightpaths            
            error (['Non Valid Wavelength in fibre ', int2str(m)]);
        end
        % Checking the Wavelength Clashing Constraint in the current fibre
        idsOfLPsInThisFibre= find(netState.lightpathRoutingMatrix(:,m)>0);
        if (numel(unique(netState.lightpathRoutingMatrix(idsOfLPsInThisFibre,m))) ~= numel(idsOfLPsInThisFibre)), % cambio (ultimo numel)
            error (['Wavelength Clashing Constraint violated in fibre ', int2str(m)]);
        end
    end
end        
    
%%%% 3. Checking the LP routing over the fibre links %%%%
if (flag_checkLpRoutingMatrix ~= 0) 
    for lpID = 1:nrOfLightpaths,    
        lpLinkIds=find(netState.lightpathRoutingMatrix (lpID,:)~=0);    
        nonTraversedYetLinksIds = lpLinkIds;
        originNode = netState.lightpathTable (lpID,2);
        destinationNode = netState.lightpathTable (lpID,3);

        if (originNode == destinationNode), error (['Lp starts and ends at the same node: ', int2str(lpID)]); end    

        currentNode = originNode;    
        while (currentNode ~= destinationNode)
            outgoingLinks = phys.outgoingLinksFromNode {currentNode};
            nextHopLinkId = intersect (outgoingLinks , nonTraversedYetLinksIds);
            if isempty(nextHopLinkId) 
                error (['Not connected Lp:  ', int2str(lpID)]);
            elseif length(nextHopLinkId)>1 
                error (['Loop in Lp:  ', int2str(lpID)]);
            end
            currentNode = phys.linkTable (nextHopLinkId,2);
            nonTraversedYetLinksIds = setdiff (nonTraversedYetLinksIds , nextHopLinkId);
        end
        if numel(nonTraversedYetLinksIds) > 0 && flag_checkLoops == 1
            if flag_removeLoops == 0 , error (['There exists no contiguous routing (maybe a cycle) in the lightpath: ' num2str(lpID)]); end
            newNetState.lightpathRoutingMatrix (lpID,nonTraversedYetLinksIds) = 0;
            disp (['Possible cycle eliminated in lp:' num2str(lpID)]);
            flag_modifiedNetState = 1;
        end        
    end
end

%%%% 4. Checking the returned number of transmitters, receivers and converters in each node  %%%%
for node=1:nrOfNodes
    outgoingLpsAttendingToNetState = find (netState.lightpathTable(:,2) == node);
    incomingLpsAttendingToNetState = find (netState.lightpathTable(:,3) == node);
    numberOutgoingLpsAttendingToNetState = length (outgoingLpsAttendingToNetState);
    numberIncomingLpsAttendingToNetState = length (incomingLpsAttendingToNetState);
    candidateTraversingLps = setdiff (1:nrOfLightpaths , union (outgoingLpsAttendingToNetState,incomingLpsAttendingToNetState));
    numberOfTraversingLpsWhichChangeLambda = 0;
    
    % this checks the wavelength converters
    if (flag_checkLpRoutingMatrix ~= 0) 
        for lpId=candidateTraversingLps
            linksTraversed = find (netState.lightpathRoutingMatrix (lpId,:) > 0);
            inputLink = intersect (phys.incomingLinksToNode {node} , linksTraversed);
            outputLink = intersect (phys.outgoingLinksFromNode {node} , linksTraversed);
            if (length (inputLink) == 1) && (length (outputLink) == 1)
                if netState.lightpathRoutingMatrix (lpId,inputLink) ~= netState.lightpathRoutingMatrix (lpId,outputLink)
                    numberOfTraversingLpsWhichChangeLambda = numberOfTraversingLpsWhichChangeLambda + 1;
                end
            end
        end
    end 
    if (netState.numberOfOccupiedTxs(node) ~= numberOutgoingLpsAttendingToNetState)
        error (['Node ' int2str(node) ': error in num TX in NetState: differ from lpTable-Routing']);
    end
    if (netState.numberOfOccupiedRxs(node) ~= numberIncomingLpsAttendingToNetState)
        error (['Node ' int2str(node) ': error in num RX in NetState: differ from lpTable-Routing']);
    end
    if (netState.numberOfOccupiedTWCs(node) ~= numberOfTraversingLpsWhichChangeLambda)
        error (['Node ' int2str(node) ': error in num TX in NetState: differ from lpTable-Routing']);
    end
    if (netState.numberOfOccupiedTxs(node) > makeRowVector(phys.numberTxPerNode(node)))
        error (['Node ' int2str(node) ': error in num TX in NetState: higher than available']);
    end
    if (netState.numberOfOccupiedRxs(node) > makeRowVector(phys.numberRxPerNode(node)))
        error (['Node ' int2str(node) ': error in num RX in NetState: higher than available']);
    end
    if (netState.numberOfOccupiedTWCs(node) > makeRowVector(phys.numberTWCPerNode(node)))
        error (['Node ' int2str(node) ': error in num TWC in NetState: higher than available']);
    end    
end
    
%%%% 5. Checking the flow routing over the lightpaths %%%%

if (min (min (netState.flowRoutingMatrix)) < 0) || (min (netState.flowTable(:,4)) < 0)
    error ('Non Positive Traffic in a flow ');
end

if not(isempty(trafficMatrix)),
    computedTrafficMatrixFromNetState = zeros(nrOfNodes);
    for flowID = 1:nrOfFlows,
        ingressNode = netState.flowTable(flowID,2);
        egressNode = netState.flowTable(flowID,3);
        computedTrafficMatrixFromNetState(ingressNode,egressNode) = computedTrafficMatrixFromNetState(ingressNode,egressNode) + netState.flowTable(flowID,4);
    end
    if ((flag_lossesAllowed)&((computedTrafficMatrixFromNetState - trafficMatrix) > 1E-6)), error('The offered traffic matrix computed from "netState.flowTable" exceeds the offered traffic from "trafficMatrix"');end
    if (not(flag_lossesAllowed)&(abs(computedTrafficMatrixFromNetState - trafficMatrix) > 1E-6)), error('The offered traffic matrix computed from "netState.flowTable" does not match the offered traffic from "trafficMatrix"');end
end

for flowID = 1:nrOfFlows,
    ingressNode = netState.flowTable(flowID,2);
    egressNode = netState.flowTable(flowID,3);
    if not(isempty(trafficMatrix)),% offered traffic computed from the "trafficMatrix"
        offeredTraffic = trafficMatrix(ingressNode,egressNode);
    else % offered traffic computed from the "netState.flowTable"
        offeredTraffic = netState.flowTable(flowID,4);
    end       
    outgoingLPsFromIngressNode = makeRowVector(find(netState.lightpathTable(:,2)==ingressNode));
    incomingLPsToEgressNode = makeRowVector(find(netState.lightpathTable(:,3)==egressNode));
    trafficOutgoingIngressNode = sum(netState.flowRoutingMatrix(flowID, outgoingLPsFromIngressNode));
    trafficIncomingEgressNode = sum(netState.flowRoutingMatrix(flowID, incomingLPsToEgressNode));
    
    if (ingressNode == egressNode), error (['Flow starts and ends at the same node: ', int2str(lpID)]); end  
    if (abs(trafficOutgoingIngressNode - trafficIncomingEgressNode) > 1E-6), 
        error (['The traffic injected in the ingress does not match the traffic received at the egress, for flow ', num2str(flowID) ]); 
    end     
    if flag_lossesAllowed,
        if (trafficOutgoingIngressNode - offeredTraffic > 1E-6), 
            error (['The traffic injected in the ingress in flow ', num2str(flowID), ' exceeds the offered traffic' ]); 
        end 
        if (trafficIncomingEgressNode - offeredTraffic > 1E-6), 
            error (['The traffic received by the egress in flow ', num2str(flowID), ' exceeds the offered traffic' ]); 
        end
    else
        if (abs(trafficOutgoingIngressNode - offeredTraffic) > 1E-6), 
            error (['The traffic injected in the ingress in flow ', num2str(flowID), ' does not match the offered traffic' ]); 
        end
        if (abs(trafficIncomingEgressNode - offeredTraffic) > 1E-6), 
            error (['The traffic received by the egress in flow ', num2str(flowID), ' does not match the offered traffic' ]); 
        end
    end
    
    traversedLpIds = find(netState.flowRoutingMatrix (flowID,:) > 0);
    traversedLpInfo = [traversedLpIds' netState.lightpathTable(traversedLpIds,2:3)];
   
    % check cycles at the ingress and at the egress
%     if (numel (find(traversedLpInfo(:,2) == egressNode) > 0)), error ('The flow traverses a lightpath outgoing from the egress node: a cycle'); end
%     if (numel (find(traversedLpInfo(:,3) == ingressNode) > 0)), error ('The flow traverses a lightpath incoming to the ingress node: a cycle'); end

    linkTable = [traversedLpInfo ones(numel(traversedLpIds),1)];
    remainingRoutedTrafficPerTraversedLp = netState.flowRoutingMatrix (flowID,traversedLpIds);
    while (1)        
        % establish shortest paths in the topology composed by the
        % traversed lps. calculate the amount of flow it can carry
        % then remove it, and repeat. at the end, all the flow should be
        % carried
        [sequenceOfSPFLinkIds,sequenceOfSPFNodeIds , totalCost] = libraryGraph_shortestPath (linkTable(:,2:end) , ingressNode , egressNode);
        if (numel(sequenceOfSPFLinkIds) == 0), break; end
        trafficRouted = min(remainingRoutedTrafficPerTraversedLp(sequenceOfSPFLinkIds));
        remainingRoutedTrafficPerTraversedLp(sequenceOfSPFLinkIds) = remainingRoutedTrafficPerTraversedLp(sequenceOfSPFLinkIds) - trafficRouted*ones(1,numel(sequenceOfSPFLinkIds));
        lpsWithoutSpareCapacity = find (remainingRoutedTrafficPerTraversedLp < 1E-6);
        linkTable(lpsWithoutSpareCapacity,:) = []; % remove the links without spare capacity
        remainingRoutedTrafficPerTraversedLp (lpsWithoutSpareCapacity) = []; % remove the links without spare capacity
    end
    
    if (numel (remainingRoutedTrafficPerTraversedLp) > 0) || (numel (linkTable)> 0) && flag_checkLoops == 1, 
        if flag_removeLoops == 0 , error ('There exists no contiguous routing (an isolated loop in the flow), or flow is not conserved'); end
        % I eliminate the traffic associated to the cycles
        idOfTheLpsRelatedToCycles = linkTable(:,1);  
        cycleFlow = min(newNetState.flowRoutingMatrix(flowID,idOfTheLpsRelatedToCycles));
        newNetState.flowRoutingMatrix (flowID,idOfTheLpsRelatedToCycles) = newNetState.flowRoutingMatrix (flowID,idOfTheLpsRelatedToCycles) - cycleFlow; 
        flag_modifiedNetState = 1;
        disp (['Cycle eliminated in flow:' num2str(flowID) ':(' num2str(netState.flowTable(flowID,2)) '-->' num2str(netState.flowTable(flowID,3)) ')']);
    end
    
    %A Lightpath Flow (sum of all the flows routed on a lightpath) must not
    %be overcome the lightpath capacity.     
    lightpathFlows = sum(netState.flowRoutingMatrix,1);
    for lpID = 1:nrOfLightpaths,
        if (lightpathFlows(lpID) - phys.lightpathCapacity) > 1E-6, error(['The flow in the lightpath ',num2str(lpID),':(' num2str(netState.lightpathTable(lpID,2)) '-->' num2str(netState.lightpathTable(lpID,3)) ') exceeds the WDM channel capacity (', num2str(phys.lightpathCapacity),' Gbps).']); end
    end
           
end