diff --git a/+af_tools/nacaairfoil.m b/+af_tools/nacaairfoil.m
index 17c4fa38de0016dab7bf0297aa98fba87fddfdc4..098bec226ef9f364a59cc38fae8730b16a4fcb8f 100644
--- a/+af_tools/nacaairfoil.m
+++ b/+af_tools/nacaairfoil.m
@@ -1,4 +1,4 @@
-function [x, y] = nacaairfoil(varargin)
+function [coord] = nacaairfoil(varargin)
     % NACAAIRFOIL Generates the full coordinates of a NACA 4 or 5 digits airfoil.
     %   By default, the airfoil is generated with 100 points distributed following a half-cosine
     %   rule in order to be finer at the leading and trailing edges. The default output has a zero
@@ -16,13 +16,13 @@ function [x, y] = nacaairfoil(varargin)
     % -----
     %
     % Usage:
-    %   [X, Y] = NACAAIRFOIL prompts the user for the airfoil denomination, then determines its
+    %   [COORD] = NACAAIRFOIL prompts the user for the airfoil denomination, then determines its
     %   coordinates.
     %
-    %   [X, Y] = NACAAIRFOIL(DIGITS) determines the coordinates of the airfoil using the airfoil
+    %   [COORD] = NACAAIRFOIL(DIGITS) determines the coordinates of the airfoil using the airfoil
     %   denomination given by DIGITS.
     %
-    %   [X, Y] = NACAAIRFOIL(DIGITS, NPOINTS) determines the NPOINTS coordinates of the airfoil
+    %   [COORD] = NACAAIRFOIL(DIGITS, NPOINTS) determines the NPOINTS coordinates of the airfoil
     %   named after DIGITS.
     %
     %   [...] = NACAAIRFOIL(..., 'spacing', SP) uses the defined spacing rule to determine to
@@ -44,10 +44,10 @@ function [x, y] = nacaairfoil(varargin)
     %   nPoints : Number of output coordinates
     %
     % Output:
-    %   [X, Y] : Airfoil coordinates, with 0 <= X <= 1.
-    %            The output coordinates follow the most usual order. They start at the trailing
-    %            edge, along the upper surface to the leading edge and back around the lower surface
-    %            to trailing edge.
+    %   coord : Airfoil coordinates, with 0 <= X <= 1.
+    %           The output coordinates follow _Selig_ order. They start at the trailing
+    %           edge, along the upper surface to the leading edge and back around the lower surface
+    %           to trailing edge.
     %
     % Example:
     %   NACAAIRFOIL
@@ -55,6 +55,8 @@ function [x, y] = nacaairfoil(varargin)
     %   NACAAIRFOIL('24012', 120)
     %   NACAAIRFOIL('24012', 120, 'spacing', 'halfcosine')
     %   NACAAIRFOIL('24012', 120, 'zerote', true)
+    %   NACAAIRFOIL('24012', 120, 'savedat', true)
+    %   NACAAIRFOIL('24012', 120, 'outputFormat', 'Lednicer')
     %
     % See also: NACACAMBER.
     %
@@ -73,22 +75,22 @@ function [x, y] = nacaairfoil(varargin)
     import af_tools.nacacamber
     import af_tools.utils.*
 
-    OPTION_LIST = {'spacing', 'zerote', 'savedat'};
+    OPTION_LIST = {'spacing', 'zerote', 'savedat', 'outputFormat'};
 
     % Extract and validate the inputs
     [digits, nPoints, idxOpts] = parsenacainputs(OPTION_LIST, varargin{:});
-    [spacing, zerote, savedat] = parseoptioninputs(varargin{idxOpts:end});
+    [spacing, zerote, savedat, outputFormat] = parseoptioninputs(varargin{idxOpts:end});
 
     % -------------------------------------------
     % Coordinates calculation
 
     % Get camberline (we use ceil((nPoints+1)/2) to ensure the airfoil has the same number of points
     % on the upper and lower surfaces.
-    [xc, yc, gradY] = nacacamber(digits, ceil((nPoints + 1) / 2), 'spacing', spacing);
+    [coord, gradY] = nacacamber(digits, ceil((nPoints + 1) / 2), 'spacing', spacing);
 
     % Flip camberline, so the first element is the trailing edge
-    xc = fliplr(xc);
-    yc = fliplr(yc);
+    xc = fliplr(coord(:, 1));
+    yc = fliplr(coord(:, 2));
 
     % Thickness value (checks were done in nacacamber)
     t = str2double(digits(end - 1:end)) / 100;
@@ -108,27 +110,31 @@ function [x, y] = nacaairfoil(varargin)
     yt = t / 0.2 * (a0 * xc.^0.5 + a1 * xc + a2 * xc.^2 + a3 * xc.^3 + a4 * xc.^4);
 
     % Upper and lower surface coordinates
-    theta = atan(gradY);
-    xu = xc - yt .* sin(theta);
-    yu = yc + yt .* cos(theta);
-    xl = xc + yt .* sin(theta);
-    yl = yc - yt .* cos(theta);
-
-    % Full coordinates
-    x = [xu, fliplr(xl(1:end - 1))];
-    y = [yu, fliplr(yl(1:end - 1))];
+    theta = vecttocol(atan(gradY));
+    xu = vecttocol(xc - yt .* sin(theta));
+    yu = vecttocol(yc + yt .* cos(theta));
+    xl = vecttocol(xc + yt .* sin(theta));
+    yl = vecttocol(yc - yt .* cos(theta));
+
+    % Fromats according to desired output
+    switch outputFormat
+        case 'Selig'
+            [coord, array] = selig([xu, yu], [xl, yl]);
+        case 'Lednicer'
+            [coord, array] = lednicer([xu, yu], [xl, yl]);
+    end
 
     % Save coordinates to dat file
     if savedat
         header = ['NACA', digits, ' - ', spacing, ' spacing'];
-        filename = ['naca', digits, '.dat'];
-        savetodat(filename, header, x, y);
+        filename = ['naca', digits, '-', outputFormat, '.dat'];
+        savetodat(filename, header, array);
     end
 
 end
 
 % --------------------------------------------------------------------------------------------------
-function [spacing, zerote, savedat] = parseoptioninputs(varargin)
+function [spacing, zerote, savedat, outputFormat] = parseoptioninputs(varargin)
     % PARSEOPTIONINPUTS Parses the options and checks their validity
 
     import af_tools.utils.*
@@ -138,19 +144,24 @@ function [spacing, zerote, savedat] = parseoptioninputs(varargin)
     DEFAULT_ZEROTE = true;
     DEFAULT_SAVEDAT = false;
     ALLOWED_SPACING = {'linear', 'cosine', 'halfcosine'};
+    DEFAULT_FORMAT = 'Selig';
+    ALLOWED_FORMATS = {'Selig', 'Lednicer'};
 
     % Option validator
     validLogical = @(x) validateattributes(x, {'logical'}, {'scalar'}, mfilename());
     validSpacing = @(x) any(validatestring(x, ALLOWED_SPACING, mfilename()));
+    validFormat = @(x) any(validatestring(x, ALLOWED_FORMATS, mfilename()));
 
     % Parse options
     p = inputParser;
     addParameter(p, 'spacing', DEFAULT_SPACING, validSpacing);
     addParameter(p, 'zerote', DEFAULT_ZEROTE, validLogical);
     addParameter(p, 'savedat', DEFAULT_SAVEDAT, validLogical);
+    addParameter(p, 'outputFormat', DEFAULT_FORMAT, validFormat);
     parse(p, varargin{:});
     spacing = p.Results.spacing;
     zerote = p.Results.zerote;
     savedat = p.Results.savedat;
+    outputFormat = p.Results.outputFormat;
 
 end
diff --git a/+af_tools/nacacamber.m b/+af_tools/nacacamber.m
index bd6a76132c0fcfb7f6653b9b138675f2c403515d..d899f936e58291cd74ded5b0a971808702c02fa8 100644
--- a/+af_tools/nacacamber.m
+++ b/+af_tools/nacacamber.m
@@ -1,4 +1,4 @@
-function [xc, yc, gradY] = nacacamber(varargin)
+function [coord, gradY] = nacacamber(varargin)
     % NACACAMBER Generates the coordinates of the camberline for a NACA 4 or 5 digits airfoil.
     %   By default, the camberline is generated with 100 points distributed following a half-cosine
     %   rule in order to be finer at the leading and trailing edges.
@@ -53,7 +53,7 @@ function [xc, yc, gradY] = nacacamber(varargin)
     narginchk(0, 4);
 
     % Import other functions from this package
-    import af_tools.utils.parsenacainputs
+    import af_tools.utils.*
 
     OPTION_LIST = {'spacing'};
 
@@ -70,6 +70,10 @@ function [xc, yc, gradY] = nacacamber(varargin)
         [xc, yc, gradY] = camber5digits(digits, nPoints, spacing);
     end
 
+    xc = vecttocol(xc);
+    yc = vecttocol(yc);
+    coord = [xc, yc];
+
 end
 
 % --------------------------------------------------------------------------------------------------
@@ -102,7 +106,7 @@ function [xc, yc, gradY] = camber4digits(digits, nPoints, spacing)
     %   P  = Position of the max camber divided by 10
     %   XX = Maximum thickness divided by 100
 
-    import af_tools.utils.spacedvector
+    import af_tools.utils.*
 
     % Airfoil parameters
     m = str2double(digits(1)) / 100;
diff --git a/CHANGELOG.md b/CHANGELOG.md
index ff3f266b97403c5479a83bc4d4179f56496981b1..a1da4f7b8dccffa767dff218e05139692de956fa 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -9,6 +9,18 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ### Added
 
+### Changed
+
+### Deprecated
+
+### Removed
+
+### Fixed
+
+## [2.0.0] - 2022-04-24
+
+### Added
+
 - Start using MISS_HIT for style and code analysis.
 - ci: Add MISS_HIT job in pipeline
 - ci: Add pre-commit hook for miss_hit
@@ -16,8 +28,12 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 ### Changed
 
 - **BREAKING**: Renamed **uiuccleaner** in **formatairfoilcoord**
+- **BREAKING**: Changed output format for **nacacamber** and **nacaairfoil**
+
 - Feat: Add possibility to select output style for **formatairfoilcoord**
 - Feat: Add Selig and Lednicer utils for airfoil formatting
+- Feat: Add option to select output format for **nacaairfoil**
+
 - Style: Adopt style from MISS_HIT
 - Copyright notice: change main copyright holder to the University
 - Doc: Revamp `CONTRIBUTING.md` to put emphasis on MISS_HIT and conventions
@@ -25,12 +41,6 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - Refactor:  Rename lift, drag and moment coefficient variables to avoid conflict
   with `cd` command.
 
-### Deprecated
-
-### Removed
-
-### Fixed
-
 ## [1.2.0] - 2022-04-10
 
 ### Added
@@ -63,7 +73,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 - Initial release
 
-[Unreleased]: https://gitlab.uliege.be/am-dept/matlab_airfoil_toolbox/-/compare/v1.2.0...master
+[Unreleased]: https://gitlab.uliege.be/am-dept/matlab_airfoil_toolbox/-/compare/v2.0.0...master
+[2.0.0]: https://gitlab.uliege.be/am-dept/matlab_airfoil_toolbox/-/compare/v1.2.0...v2.0.0
 [1.2.0]: https://gitlab.uliege.be/am-dept/matlab_airfoil_toolbox/-/compare/v1.1.0...v1.2.0
 [1.1.0]: https://gitlab.uliege.be/am-dept/matlab_airfoil_toolbox/-/compare/v1.0.0...v1.1.0
 [1.0.0]: https://gitlab.uliege.be/am-dept/matlab_airfoil_toolbox/-/releases/v1.0.0
diff --git a/README.md b/README.md
index 7e181b15db0024bd1ecfa1852a6ab017b988e497..11ce282b24f4162fa7e062c32d2625f986d19d0d 100644
--- a/README.md
+++ b/README.md
@@ -167,9 +167,9 @@ points for each surface.
 
 ```matlab
 import af_tools.xf2mat
-Polar = uiuccleaner
-Polar = uiuccleaner(inputDir, inputFiles, 'autosave', true)
-Polar = uiuccleaner(inputDir, inputFiles, 'overwrite', true)
+Af = formatairfoilcoord
+Af = formatairfoilcoord(inputDir, inputFiles, 'autosave', true)
+Af = formatairfoilcoord(inputDir, inputFiles, 'overwrite', true)
 ```
 
 | Input          | Example                                        | Default |
@@ -198,10 +198,10 @@ The output is ordered from the leading edge to the trailing edge.
 
 ```matlab
 import af_tools.nacacamber
-[xc, yc, gradY] = nacacamber
-[xc, yc, gradY] = nacacamber(digits)
-[xc, yc, gradY] = nacacamber(digits, nPoints)
-[xc, yc, gradY] = nacacamber(digits, nPoints, 'spacing', 'halfcosine')
+[coord, gradY] = nacacamber
+[coord, gradY] = nacacamber(digits)
+[coord, gradY] = nacacamber(digits, nPoints)
+[coord, gradY] = nacacamber(digits, nPoints, 'spacing', 'halfcosine')
 ```
 
 | Options   | Example                            | Default        |
@@ -217,11 +217,11 @@ surface to trailing edge.
 
 ```matlab
 import af_tools.nacaairfoil
-[x, y] = nacaairfoil
-[x, y] = nacaairfoil(digits)
-[x, y] = nacaairfoil(digits, nPoints)
-[x, y] = nacaairfoil(digits, nPoints, 'spacing', spacing, 'zerote', true)
-[x, y] = nacaairfoil(digits, nPoints, 'savedat', true)
+[coord] = nacaairfoil
+[coord] = nacaairfoil(digits)
+[coord] = nacaairfoil(digits, nPoints)
+[coord] = nacaairfoil(digits, nPoints, 'spacing', spacing, 'zerote', true)
+[coord] = nacaairfoil(digits, nPoints, 'savedat', true)
 ```
 
 | Options   | Example                            | Default        |
diff --git a/tests/test_nacaairfoil.m b/tests/test_nacaairfoil.m
index b3fe143474b032fc5ee4597fecf7ded7c1cca2a2..6b14f2ad03028f031e230ec19ff71a0ba325c8a5 100644
--- a/tests/test_nacaairfoil.m
+++ b/tests/test_nacaairfoil.m
@@ -170,19 +170,15 @@ function test_correctNbCoordinates(testCase)
     evenNPoints = 240;
     oddNPoints = 241;
 
-    [xc4even, yc4even] = af_tools.nacaairfoil('0012', evenNPoints);
-    [xc5even, yc5even] = af_tools.nacaairfoil('24012', evenNPoints);
-    [xc4odd, yc4odd] = af_tools.nacaairfoil('0012', oddNPoints);
-    [xc5odd, yc5odd] = af_tools.nacaairfoil('24012', oddNPoints);
-
-    verifyEqual(testCase, length(xc4even), evenNPoints + 1);
-    verifyEqual(testCase, length(yc4even), evenNPoints + 1);
-    verifyEqual(testCase, length(xc5even), evenNPoints + 1);
-    verifyEqual(testCase, length(yc5even), evenNPoints + 1);
-    verifyEqual(testCase, length(xc4odd), oddNPoints);
-    verifyEqual(testCase, length(yc4odd), oddNPoints);
-    verifyEqual(testCase, length(xc5odd), oddNPoints);
-    verifyEqual(testCase, length(yc5odd), oddNPoints);
+    [coord4even] = af_tools.nacaairfoil('0012', evenNPoints);
+    [coord5even] = af_tools.nacaairfoil('24012', evenNPoints);
+    [coord4odd] = af_tools.nacaairfoil('0012', oddNPoints);
+    [coord5odd] = af_tools.nacaairfoil('24012', oddNPoints);
+
+    verifyEqual(testCase, size(coord4even, 1), evenNPoints + 1);
+    verifyEqual(testCase, size(coord5even, 1), evenNPoints + 1);
+    verifyEqual(testCase, size(coord4odd, 1), oddNPoints);
+    verifyEqual(testCase, size(coord5odd, 1), oddNPoints);
 
 end
 
@@ -197,17 +193,17 @@ function test_correctSpacing(testCase)
     betacos = linspace(0, pi / 2, ceil((nPoints + 1) / 2));
     cosine = 1 - cos(betacos);
 
-    [xc_lin, ~] = af_tools.nacaairfoil('0012', nPoints, 'spacing', 'linear');
-    [xc_cos, ~] = af_tools.nacaairfoil('0012', nPoints, 'spacing', 'cosine');
-    [xc_half, ~] = af_tools.nacaairfoil('0012', nPoints, 'spacing', 'halfcosine');
+    [coord_lin] = af_tools.nacaairfoil('0012', nPoints, 'spacing', 'linear');
+    [coord_cos] = af_tools.nacaairfoil('0012', nPoints, 'spacing', 'cosine');
+    [coord_half] = af_tools.nacaairfoil('0012', nPoints, 'spacing', 'halfcosine');
 
     % Test the two halves of the output vector (i.e. upper and lower surfaces)
-    verifyEqual(testCase, xc_lin(1:(end + 1) / 2), fliplr(linear));
-    verifyEqual(testCase, xc_cos(1:(end + 1) / 2), fliplr(cosine));
-    verifyEqual(testCase, xc_half(1:(end + 1) / 2), fliplr(halfcos));
-    verifyEqual(testCase, xc_lin(end:-1:(end + 1) / 2), fliplr(linear));
-    verifyEqual(testCase, xc_cos(end:-1:(end + 1) / 2), fliplr(cosine));
-    verifyEqual(testCase, xc_half(end:-1:(end + 1) / 2), fliplr(halfcos));
+    verifyEqual(testCase, coord_lin(1:(end + 1) / 2, 1), fliplr(linear)');
+    verifyEqual(testCase, coord_cos(1:(end + 1) / 2, 1), fliplr(cosine)');
+    verifyEqual(testCase, coord_half(1:(end + 1) / 2, 1), fliplr(halfcos)');
+    verifyEqual(testCase, coord_lin(end:-1:(end + 1) / 2, 1), fliplr(linear)');
+    verifyEqual(testCase, coord_cos(end:-1:(end + 1) / 2, 1), fliplr(cosine)');
+    verifyEqual(testCase, coord_half(end:-1:(end + 1) / 2, 1), fliplr(halfcos)');
 
 end
 
@@ -218,10 +214,9 @@ function test_correctDefaults(testCase)
     DEF_SPACING = 'halfcosine';
     DEF_ZEROTE = true;
 
-    [xc, yc] = af_tools.nacaairfoil('0012');
-    [xc_def, yc_def] = af_tools.nacaairfoil('0012', DEF_NPOINTS, 'spacing', DEF_SPACING, 'zerote', DEF_ZEROTE);
+    [coord] = af_tools.nacaairfoil('0012');
+    [coord_def] = af_tools.nacaairfoil('0012', DEF_NPOINTS, 'spacing', DEF_SPACING, 'zerote', DEF_ZEROTE);
 
-    verifyEqual(testCase, xc, xc_def);
-    verifyEqual(testCase, yc, yc_def);
+    verifyEqual(testCase, coord, coord_def);
 
 end
diff --git a/tests/test_nacacamber.m b/tests/test_nacacamber.m
index ca067974f2e462acfe72f97b4d42dc5219cedfd5..417bc60d26c6a9035031a56debaa7c6c797375a8 100644
--- a/tests/test_nacacamber.m
+++ b/tests/test_nacacamber.m
@@ -182,13 +182,13 @@ function test_correctSpacing(testCase)
     betacos = linspace(0, pi / 2, nPoints);
     cosine = 1 - cos(betacos);
 
-    [xc_lin, ~] = af_tools.nacacamber('0012', nPoints, 'spacing', 'linear');
-    [xc_cos, ~] = af_tools.nacacamber('0012', nPoints, 'spacing', 'cosine');
-    [xc_half, ~] = af_tools.nacacamber('0012', nPoints, 'spacing', 'halfcosine');
+    [coord_lin] = af_tools.nacacamber('0012', nPoints, 'spacing', 'linear');
+    [coord_cos] = af_tools.nacacamber('0012', nPoints, 'spacing', 'cosine');
+    [coord_half] = af_tools.nacacamber('0012', nPoints, 'spacing', 'halfcosine');
 
-    verifyEqual(testCase, xc_lin, linear);
-    verifyEqual(testCase, xc_cos, cosine);
-    verifyEqual(testCase, xc_half, halfcos);
+    verifyEqual(testCase, coord_lin(:, 1), linear');
+    verifyEqual(testCase, coord_cos(:, 1), cosine');
+    verifyEqual(testCase, coord_half(:, 1), halfcos');
 
 end
 
@@ -218,18 +218,18 @@ function test_correctMaxCamber(testCase)
     maxCambPos5 = 0.2;
 
     % Calculate camberline
-    [xc4, yc4] = af_tools.nacacamber(naca4, 100000, 'spacing', 'linear');
-    [xc5, yc5] = af_tools.nacacamber(naca5, 100000, 'spacing', 'linear');
+    [coord4] = af_tools.nacacamber(naca4, 100000, 'spacing', 'linear');
+    [coord5] = af_tools.nacacamber(naca5, 100000, 'spacing', 'linear');
 
     % Get max camber value and position from camberline coordinates
-    [maxVal4, maxPos4] = max(yc4);
-    [~, maxPos5] = max(yc5);
+    [maxVal4, maxPos4] = max(coord4(:, 2));
+    [~, maxPos5] = max(coord5(:, 2));
 
     % Naca 4
     verifyEqual(testCase, maxVal4, maxCambVal4, 'RelTol', 1e-3);
-    verifyEqual(testCase, xc4(maxPos4), maxCambPos4, 'RelTol', 1e-3);
+    verifyEqual(testCase, coord4(maxPos4, 1), maxCambPos4, 'RelTol', 1e-3);
 
     % Naca 5
-    verifyEqual(testCase, xc5(maxPos5), maxCambPos5, 'RelTol', 1e-3);
+    verifyEqual(testCase, coord5(maxPos5, 1), maxCambPos5, 'RelTol', 1e-3);
 
 end