diff --git a/read.m b/read.m
index f18cf222c35c71b90bb3f46aedec933a64d2aa73..377652b5fa19cf25944549dcec83e51e64873d71 100644
--- a/read.m
+++ b/read.m
@@ -1,19 +1,25 @@
 function frames = read(varargin)
-%READ Summary of this function goes here
-%   Detailed explanation goes here
-
+%READ Write yuv file from a cell array of frames
+%   READ(frames,'name',name) read yuv file with corresponding to name filename
+%   READ(frames,'recDir',recDir) read from directory recDir
+%   READ(frames,'subSamp',subSamp) specify YUV sampling scheme
+%   READ(frames,'writePrecision',writePrecision) specify write data type
+%   READ(frames,'readPrecision',readPrecision) specify read data type
+
+%% Parse input parameters
 subSamps  = {'400','420','422','444'};
 types = {'bit','ubit','int','uint','float','single','double'};
 validatePrecision = @(prec) any(validatestring(utils.precision(prec),types));
 defRecDir = fullfile(pwd,'Rec');
 
 p = inputParser; p.StructExpand = true; p.KeepUnmatched = true;
-p.addParameter('recDir', defRecDir , @ischar);
 p.addParameter('name'  , 'sequence', @ischar);
+p.addParameter('recDir', defRecDir , @ischar);
 
 p.parse(varargin{:});
 yuvParams = p.Results;
 
+%% Set remaining default parameters from matfile if it exists
 try
     m = get_mat(yuvParams.recDir,yuvParams.name,false);
     defSubSamp = m.subSamp;
@@ -26,27 +32,26 @@ catch
 end
 
 p.addParameter('subSamp'       , defSubSamp       , @(x) any(validatestring(x,subSamps) ));
-p.addParameter('readPrecision' , defReadPrecision , validatePrecision);
 p.addParameter('writePrecision', defWritePrecision, validatePrecision);
+p.addParameter('readPrecision' , defReadPrecision , validatePrecision);
 
 p.parse(varargin{:});
 yuvParams = p.Results;
 
-subSamp = yuvParams.subSamp;
-readPrecision = yuvParams.readPrecision;
-writePrecision = yuvParams.writePrecision;
 name = yuvParams.name;
 recDir  = yuvParams.recDir;
+subSamp = yuvParams.subSamp;
+writePrecision = yuvParams.writePrecision;
+readPrecision = yuvParams.readPrecision;
+precision = [writePrecision '=>' readPrecision];
 
+%% Recover yuv size
 m = get_mat(recDir,name,false);
 yuvSize = m.yuvSize;
 imgRes  = yuvSize(1:2);
 imgSize = yuvSize(3:end);
 
-precision = [writePrecision '=>' readPrecision];
-
-filename = fullfile(recDir,[name,'.yuv']);
-
+%% YUV 400/420/422: define interpolation sampling grids
 [xq,xgv] = deal(1:imgSize(1));
 [yq,ygv] = deal(1:imgSize(2));
 
@@ -62,22 +67,30 @@ subImgSize = [numel(xgv),numel(ygv)];
 
 imgGrid = zeros(imgRes);
 
+%% Read frames iteratively
+filename = fullfile(recDir,[name,'.yuv']);
 fileID = fopen(filename,'r');
+
 if strcmp(subSamp,'400')
     frames = arrayfun(@readY,imgGrid,'UniformOutput',false);
 else
     frames = arrayfun(@readYUV,imgGrid,'UniformOutput',false);
 end
+
 fclose(fileID);
 
+%% Auxiliary function
+    % Read a single frame of given dimensions
     function frame = readFrame(imgSize)
         frame = fread(fileID,flip(imgSize),precision)';
     end
 
+    % Read a single Y frame
     function Y = readY(~)
         Y = readFrame(imgSize);
     end
 
+    %Read Y, U and V frames
     function YUV = readYUV(~)
         Y = readFrame(imgSize);
         U.Values = readFrame(subImgSize);
diff --git a/write.m b/write.m
index 495f33610e56946ac0836f4e20c5c0ac62b218ad..c8d2f449a388f1d22221d0472e013528c5b6ae61 100644
--- a/write.m
+++ b/write.m
@@ -1,64 +1,74 @@
 function yuvSize = write(frames,varargin)
-%WRITE Summary of this function goes here
-%   Detailed explanation goes here
-
+%WRITE Write yuv file from a cell array of frames
+%   WRITE(frames,'name',name) use name as yuv filename
+%   WRITE(frames,'refDir',refDir) write in directory refDir
+%   WRITE(frames,'subSamp',subSamp) specify YUV sampling scheme
+%   WRITE(frames,'writePrecision',writePrecision) specify write data type
+%   WRITE(frames,'readPrecision',readPrecision) specify read data type
+
+%% Parse input parameters
 subSamps  = {'400','420','422','444'};
 types = {'bit','ubit','int','uint','float','single','double'};
 validatePrecision = @(prec) any(validatestring(utils.precision(prec),types));
 defRefDir = fullfile(pwd,'Ref');
 
 p = inputParser; p.StructExpand = true; p.KeepUnmatched = true;
-p.addParameter('refDir', defRefDir , @ischar);
 p.addParameter('name'  , 'sequence', @ischar);
+p.addParameter('refDir', defRefDir , @ischar);
 
 p.parse(varargin{:});
 yuvParams = p.Results;
 
+%% Set remaining default parameters from matfile if it exists
 try
     m = get_mat(yuvParams.refDir,yuvParams.name,false);
     defSubSamp = m.subSamp;
-    defReadPrecision = m.readPrecision;
     defWritePrecision = m.writePrecision;
+    defReadPrecision = m.readPrecision;
 catch
     defSubSamp = subSamps{1};
-    defReadPrecision = class(frames{1});
     defWritePrecision = class(frames{1});
+    defReadPrecision = class(frames{1});
 end
 
 p.addParameter('subSamp'       , defSubSamp       , @(x) any(validatestring(x,subSamps) ));
-p.addParameter('readPrecision' , defReadPrecision , validatePrecision);
 p.addParameter('writePrecision', defWritePrecision, validatePrecision);
+p.addParameter('readPrecision' , defReadPrecision , validatePrecision);
 
 p.parse(varargin{:});
 yuvParams = p.Results;
 
+name = yuvParams.name;
+refDir  = yuvParams.refDir;
 subSamp = yuvParams.subSamp;
 readPrecision = yuvParams.readPrecision;
 writePrecision = yuvParams.writePrecision;
-name = yuvParams.name;
-refDir  = yuvParams.refDir;
 
 imgRes  = size(frames);
 imgSize = size(frames{1});
 yuvSize = [imgRes,imgSize];
 
+%% Save parameters to read yuv file back
 m = get_mat(refDir,name,true);
-m.yuvSize = yuvSize;
-m.subSamp = subSamp;
-m.readPrecision = readPrecision;
-m.writePrecision = writePrecision;
 m.name = name;
 m.refDir  = refDir;
+m.subSamp = subSamp;
+m.writePrecision = writePrecision;
+m.readPrecision = readPrecision;
+m.yuvSize = yuvSize;
 
+%% Incompatibility (needs better handling)
 if (~strcmp(subSamp,'400')&&numel(imgSize)~=3)
     error('Incompatible frame size and subsampling');
 end
 
+%% YUV 400: pick first channel for every frame
 frames = cellfun(@double,frames,'UniformOutput',false);
 if strcmp(subSamp,'400')
     frames = cellfun(@(f) f(:,:,1),frames,'UniformOutput',false);
 end
 
+%% YUV 420/422: define interpolation sampling grids
 [xq,xgv] = deal(1:imgSize(1));
 [yq,ygv] = deal(1:imgSize(2));
 
@@ -71,31 +81,38 @@ end
 
 [U,V] = deal(griddedInterpolant({xgv,ygv},zeros(imgSize(1:2))));
 
-if ~exist(refDir,'dir')
-    mkdir(refDir);
-end
-
+%% Write frames iteratively
 filename = fullfile(refDir,[name,'.yuv']);
+
+% Checker wether yuv file already exists
 if ~exist(filename,'file')
-    [pathstr,~,~] = fileparts(filename);
-    if ~exist(pathstr,'dir')
-        mkdir(pathstr)
-    end
-    fileID = fopen(filename,'w+');
-    if strcmp(subSamp,'400')
-        cellfun(@writeFrame ,frames,'UniformOutput',false);
-    else
-        cellfun(@writeFrames,frames,'UniformOutput',false);
+        mkdir(refDir);
+else
+    if isunix && ~ismac 
+        warning('No overwriting allowed by default on Linux');
+        warning('This is to avoid conflict between concurrent threads on IGRIDA');
+        warning('You can modify this behaviour by commenting this condition in the source code');
+        return;
     end
-    fclose(fileID);
+end
+
+fileID = fopen(filename,'w+');
+
+if strcmp(subSamp,'400')
+    cellfun(@writeFrame ,frames,'UniformOutput',false);
 else
-    warning('YUV file already exists');
+    cellfun(@writeFrames,frames,'UniformOutput',false);
 end
 
+fclose(fileID);
+
+%% Auxiliary functions
+    % Write a single Y frame assuming 400 sampling scheme
     function writeFrame(frame)
         fwrite(fileID,cast(frame',readPrecision),writePrecision);
     end
 
+    % Write Y,U and V frames assuming 422/420 sampling scheme
     function writeFrames(frames)
         writeFrame(frames(:,:,1));
         U.Values = frames(:,:,2);