A. Reading data from a file for plotting and exploration : Given a (plain text!) file containing N columns of values, we need to read the values into matlab for exploration and plotting. To do this, we need to know the name of the file, and how data is arranged in it, so need to view it in (a plain text!) editor. 1. Write a script "readFile.m" which: a. reads a given file "file.dat" into a matrix "mat", b. extracts the columns of the data from "mat" into arrays "c1", "c2", ..., "cN", c. plots some column versus another column, with a title, and labeled axes. To test and debug your code, create a "file.dat" with 3 rows and N=2 columns ( so a 32 matrix ). The first column should be sorted. You can just enter some numbers in "file.dat" in an editor. Next, test it on a 43 array to make sure it also works. Insert brief comments in you code to document what is being done. 2. To apply this on real data, I am providing a data file, produced by a code that implements a mathematical model for bacterial chemotaxis (by numerical solution of a system of Partial Differential Equations). The file prof300s.dat contains profiles (values at a grid of locations x(i) at some particular time, here at 300sec) of several quantities, labeled by column at the top ( x(i) C(i) B(i) etc ). The last column is the index "i" of x(i). Of main interest is the 3rd column, labeled "B(i)". It represents bacterial density at location x(i), at the output time (normalized, B(i) means B(i)10⁹ cells/mL). It has a "bump", and want to find its FWHM, i.e. the width of the interval [xLeft, xRight] where a horizontal line at half-maximum height intersects the bump. Biologically, this is the width of the "band" the bacteria form during chemotaxis. B. Plotting 1. Download the file prof300s.dat and save it as "prof300s.dat" . View it in an editor to see its contents. The first two lines contain labels, need to be read over with fgetl . The rest of the lines are numerical values to be scanned into a matrix, say "mat". You can see what the B-profile looks like easily in gnuplot: gnuplot> plot "prof300s.dat" u 1:3 w lines lw 2 2. Copy your "readFile.m" to "Lab8.m" and modify appropriately to achieve the following: a. Read the data from "prof300s.dat" into a matrix "mat". b. Extract each column (which is now a row of "mat") into a vector (with sensible names: x , C , B , ... ). Only need x and B, so may extract only those. Suppress output with ";", they are long... c. Find Bmax = maximum value of B, and at which index it occurs, iBmax. Set cut = Bmax/2 . Print Bmax and cut nicely (with fprintf). Run your code to make sure it produces values of Bmax, cut. d. Plot B vs x, as a red curve. Label the axes. The units of x are micrometers (m). B is dimensionless. Run your code to make sure it produces the plot. e. To see where we are cutting the bump, also plot the line at height cut, on the same plot, as a black horizontal line. [Note: Once you know the value of "cut", say 1.5, it's easy to do it in Gnuplot (outside of Matlab): gnuplot> cut=1.5; plot "prof300s.dat" u 1:3 w lines lw 2 , cut lw 3 ] To do this in Matlab, will need to create an array, call it "cutArray", of same length as x, with constant values = cut. Again, run your code to make sure it produces the plot with both the bump and the cut-line. f. Is this value of cut a good choice here? does it really represent half-maximum of the bump ? (note that the "base" value of B is 1, not zero). Figure out a better expression for setting "cut". C. Find FWHM numerically, from arrays of points. 1. Create a Matlab function in a file "fwhm.m", which accepts two arrays (x, y) as input, and returns FWHM of y, by implementing the following algorithm (discussed in Part 5 of HW1).

Find the maximum value Peak of y, and the index where it occurs, iPeak .

Set the cutoff level at half maximum: cut = Peak/2, or more appropriately from B.2.f.

Find the first point x_i (before iPeak) at which y_i cut. Call it xLeft.

Find the first x_i (after iPeak) at which y_i cut. Call it xRight.

Then FWHM = xRight xLeft . 2. Modify your "Lab8.m" code to: a. Call the fwhm function with arguments (x, B) to find the FWHM of the B array. b. Then print it out (nicely, with fprintf) to say: "FWHM of bacteria at time 300 sec is" FWHM_value "microns, at height" cut_value