Tuesday, November 27, 2012

Exercise 1-12 Write a program that prints its input one word per line

Kernighan & Ritchie - The  C  Programming Language (Second Edition)

Chapter 1

Exercise 1-12 Write a program that prints its input one word per line

This exercise is very interesting, because it gives you an important piece of information about how the input and output are read and stored.

The solution:

#include <stdio.h>
int main () {   

    int c;
    while((c = getchar()) != EOF)
        if(c == ' ' || c == '\t' || c == '\n')
        // flush the output!!!
        // store the output!!!

return 0;

So how this program really works?
In order to figure out, you have to know something that is NOT reported in K&R2.

Every time you type something into the console, this IS NOT read immediately by your program.
Instead, it is stored in a buffer memory until a special input is inserted!
This input is the newline.

The newline tells the Operating System (NOT YOUR PROGRAM) something like: "OK dude I typed everything I wanted and NOW pass it to my program so it can process it."

The same is also true for the output. Nothing really happens on the console until the special character '\n' is given to the operating system.

If you have understood my highly informal explanation, now you should be able to figure out how my solution works.

You type all your text, including words, tabs, blanks and, eventually, you press ENTER on your keyboard.

When the Operating System receives the newline command, it sends everything to your program that is now able to process, character by character, all your input.

During the processing, every character that is NOT a blank, a tab or a newline, is stored in a buffer for future release to the console. On the other hand, when the program encounters a blank, a newlines or a tab it gives the command ENTER to the operating system, that flushes ALL the content that has been stored previously in the buffer.

That's it.

For more information, please refer to A. Koenig and B.E. Moo - Accelerated C++ Practical Programming be Example - Chapter 1.

As always, if you like this post .... consider to offer me a coffee ;)

Exercise 1-7. Write a program to print the value of EOF

Kernighan Ritchie - The C Programming Language (Second Edition)
Esercise 1-7
Write a program to print the value of EOF

#include <stdio.h>
int main () {
    printf("\nValue of EOF: %d\n", EOF);
    return 0;

The result in my case is -1 .

To insert an End Of File (EOF) on a Mac console you need to press "CTRL + D"

Friday, November 9, 2012

How to fit IMPLICIT functions with Matlab.

In Matlab you can perform implicit curve fit using either "nlinfit" or "lsqcurvefit" (the last one only if you have the optimization toolbox) but you need to re-write the implicit function into an explicit one solving it numerically (i.e. using "fzero").
THERE IS NO WAY to fit directly using an implicit expression!

To learn how to successfully fit using an implicit function you can follow the instructions on this Matlab page page or to read and use my worked-out example.

I will show you how my implementation of the Langevin function with a field Weiss (that is an implicit equation). This example is interesting by itself becasue it addresses the topic of how to model an hysteresis loop.
The equation is taken from the paper by D.C. Jiles and D.L. AthertonFerromagnetic Hysteresis, IEEE Transaction on Magnetics, Vol. Mag-19, No. 5, September 1983.

For a general reference about the Langevin equation (and generally for magnetism) I recommend the nice (and inexpensive!) book by S. BlundellMagnetism in Condensed Matter.

So let's start!

So we want to create a fitting function in Matlab using this IMPLICIT equation:

As you see the above equation is implicit due to the term alpha*M.

The first step is to transform the equation into an EXPLICIT one. As I wrote in the beginning, it is necessary just to use the fzero function.
I named this function  LangevinWithWeissField.m (if you are lazy, you can download the file .m directly from this link)

function y=LangevinWithWeissField(p, x)% Implicit function: Langevin function with a Weiss field.
% y = A*(coth(B*x(i)+C*y) - 1./(B*x(i)+C*y))
% p is the parameter vector
% Physical meaning of the parameters
% x is the magnetic field H
% y is the magnetization M corresponding to n*<m>z where n is the number of
% magnetic moments,  and <m>z is the projected magnetic moment along the
% field direction (called z) of the total moment mu
%A = amplitude
%B = mu0 m/(kB*T) where kB is the Boltzman constant and T is the temperature
%C is a factor giving the Weiss field (that is proportional to the total
%magnetization y, i.e. M=n*<m>z, note that this term includes also the previous factor
% mu0 m /(kB*T)

% assign the parameters ...


y=zeros(size(x)); % define a vector to allocate the magnetization values
NN=length(x); % total length of the field vector x, i.e. B

opt = optimset('display','off');
% I out off all the messages coming from fzero. If something goes wrong, change this option to 'off'
% to see at which x values fzero failed.

for i=1:NN   
   y(i)=fsolve(@(y) y - A*(coth(B*x(i)+C*y) - 1./(B*x(i)+C*y)) , 0.0001, opt);
   % Here 0.0001 is our starting point to find the solution around 0.

end % close the function

Once you have saved the code in your Matlab working directory, all the work it's done!

You can plot the Langevin function with the parameter you like, for example (from the Matlab command window)

xx=linspace(-30, 30, 200); % magnetic field
p=[10, 2, 0.1]; %parameters
yy=LangevinWithWeissField(p, xx); % magnetization values
plot(xx, yy, '*-r'); % plot

Now let's see how to "fit" the y-data we just created, i.e. the vector yy.
I use the  LSQCURVEFIT function as done in this mathworks' link.
The code to use in the command line is

lsqcurvefit(@(params, xdata) LangevinWithWeissField(params, xdata),[1 1 1], xx, yy)

In the previous line, the vector  [1 1 1] represents the initial guess of the fitting parameters.
After a while (on my computer it took 54 seconds) the lsqcurvefit gives you the optimal fit parameters, that corresponds to the original values stored in the p vector, namely [10 2 0.1].

To conclude, here is two things to have to keep in mind.
1) It's VERY important to use a good set of values as initial guess for your fit.  With bad (=random) combination of values, your fit could not converge!

2) Using implicit functions can very easily make the fitting SLOW, terribly SLOW!
So, try to avoid their use whenever it is possible.

Any question or comment, as always, is welcome!

You could also be interested in these other posts

(A)  How to hold fitting parameters in nlinfit (MATLAB) 
(B)  edufit: a Matlab data fitting interface

I hope that you liked this post and if you did, consider to support my blog! ;) 

Thursday, November 8, 2012

How to hold parameters in nlinfit (Matlab)

If you are using nlinfit to fit your data in Matlab, maybe you are interested (as I was) in finding an easy way to hold some of your parameters (by default there no option to do it).
So how to hold, how to fix some of the parameters?

Looking on the Internet I found this "official" solution link coming directly from the Mathworks support.
It consists in re-writing a new version of nlinfit, called nlinfitsome, that fits only some of parameters. 
The solution is pretty nice and clever but it has a drawback.

The reason I use nlinfit over lsqcurvefit is because it makes relatively easy to get the error bars, however with the proposed solution this advantage is lost so I modified nlinfitsome to address this specific problem.

My function, called mod_nlinfit, has exactly the same arguments of nlinfit and gives you as output the parameters and their error bars. I know that nlinfit has a richer output, but for my needs, having the parameters and error bars is enough.

The code is quite straightforward.

The hold parameters will have no uncertainty.

Mathworks support

Matlab code:

function [p, r,  err]=mod_nlinfit(fixed, x, y, fun, beta0, varargin)
% [p, err]=mod_nlinfit(fixed, x, y, fun, beta0, varargin)
% p: vector containing the fitting results

% r: residuals
% err: error bars of the parameters
% fixed is a vector of the same length of beta0 and containing logical elements (i.e. true and
% false). True values mean FIXED parameters.
% If you have three parameters, a vector fixed [true true false]
% means that only the third parameter will change

% all the values in beta0 are copied in p however only the free (NOT fixed)
% ones will change


% fixed parameters will get zero uncertainity


% extraction of the free parameters


[p1, r1, JJ1] = nlinfit(x,y,@nested_fun,beta1,varargin{:});

% getting the error bars for the free parameters
delta=nlparci(p1, r1, JJ1);
err1=(delta(:, 2)- delta(:, 1))/2;

% update the parameters


% ... and the error bars

    % Nested function takes just the parameters to be estimated as inputs
    % It inherits the following from the outer function:
    %   fixed = logical index for fixed elements
    %   beta0 = original guess for the parameters

    function yy=nested_fun(beta, x)
        yy=fun(b, x);


This code is used in edufit

Thursday, October 11, 2012

How to check if the input is a number (Matlab)

If you are writing a code where the user is asked to give a number, it can be a good idea to check whether or not the inserted value was really a number instead of something else like a letter or some strange character (%, / _ , etc).

The easiest way to do it is asking for a string input and then attempting (via str2num) to make the conversion string to number. If the conversion fails, then the input wasn't a number!

Here a very simple code that does the job.


function [varargout]=CheckIfNumber(s)
% s is a string
% varargout is the converted number or row vector


if (isempty(num)==0)




Wednesday, October 3, 2012

How to get a timestamp (date + time) in Matlab

To complete my  Matlab interface for data fitting, I  needed a Timestamp function, namely a function combining date and time into a single string.

This can be useful when you want to save files having the same base "name" (e.g. foo-x11.dat, foo-x12.dat,  etc) but created at different times.

In Matlab the timestamp is obtained with the command

>>datestr(clock, 0)

which output is (at the time of this post!!!)

10-nov-2012 14:10:20

Unfortunately this string cannot be used into a filename due to the presence of colons.
In addition I would like to have a string starting with "year-month-day" because it is a much more convenient  way to name files created at different times.

My version of TimeStamp that produces a string in the format :


as for example


Here is the code! Enjoy!
And keep following this blog and my other one.  It's totally worth it!


function [s]=TimeStamp
% time stamp in the format
% year-month-day-hours-h-minutes-m-second-s

% Get current time as date vector


% convert date vector into the string: year-month-day
yyyy_mm_dd=datestr(t, 29);

% convert date vector into the string: day-month-year[blank]hours:minutes:seconds
s=datestr(t, 0);

% get index of the blank
index=find(s==' ');

% remove the string that is before the blank

% get the indices of the colons ":"

%change the first colon into "h" and the second one into "m". Add an "s" at the end of the

h_m_s=strcat(  s(1:(ii(1)-1)), 'h', s( (ii(1)+1): (ii(2)-1)), 'm', s( (ii(2)+1):end), 's');

% concate the two strings
s=strcat(yyyy_mm_dd, '-', h_m_s);

Alternatively, you can use the very short code provided by the first commenter :)
In this case:

function s=TimeStamp

Or the code built on the suggestion provided by the another commenter :)

function s=TimeStamp

tmp = strrep(datestr(clock), ':' , '-' );
s= strrep(tmp, ' ' , '-' );
s=strcat(s, 's');

Tuesday, October 2, 2012

Make Matlab listing ONLY subfolders cointained into a folder

From the Matlab Command Window when you type "dir" you get a list of files and folders contained in your working directory.
The question is: how can you get ONLY the folders?

Here is the solution from stackoverflow.

In the following function the output "nameFolds" is a cell variable type!!!


function nameFolds=ListSubfolders(pathFolder)

d = dir(pathFolder);
isub = [d(:).isdir];
nameFolds = {d(isub).name}';
nameFolds(ismember(nameFolds,{'.','..'})) = [];