# Using ReadList and Skip to read a file with periodic alternation of useful and useless data

I am trying a read a file that's arranged such that there's alternating set of rows containing data that I don't need and data that I need. So for my case the first n lines in the files are useless to me while the next m lines are data I need, and the next n lines are useless etc. Here is a sample.

TER
ENDMDL
REMARK    GENERATED BY TRJCONV
TITLE     Neat THF t=   6.00000
REMARK    THIS IS A SIMULATION BOX
CRYST1   32.503   32.503   32.503  90.00  90.00  90.00 P 1           1
MODEL       16
ATOM      1  O1  THF     1      19.380  32.070  27.250  1.00  0.00
ATOM      2  CA1 THF     1      19.090  30.880  26.540  1.00  0.00
ATOM      3  CA2 THF     1      19.460  33.170  26.380  1.00  0.00
ATOM      4  CB1 THF     1      18.940  31.200  25.060  1.00  0.00
ATOM      5  CB2 THF     1      19.190  32.710  24.950  1.00  0.00
ATOM      6  T21 THF     1      18.950  30.900  25.460  1.00  0.00
ATOM      7  T22 THF     1      19.290  33.050  25.300  1.00  0.00
ATOM      8  O1  THF     2      31.020  24.150   7.910  1.00  0.00
ATOM      9  CA1 THF     2      31.310  24.450   9.250  1.00  0.00
ATOM     10  CA2 THF     2      31.450  22.850   7.570  1.00  0.00
ATOM     11  CB1 THF     2      32.010  23.250   9.900  1.00  0.00
ATOM     12  CB2 THF     2      32.100  22.200   8.790  1.00  0.00
ATOM     13  T21 THF     2      31.810  23.710   9.890  1.00  0.00
ATOM     14  T22 THF     2      31.950  22.220   8.310  1.00  0.00
TER
ENDMDL
REMARK    GENERATED BY TRJCONV
TITLE     Neat THF t=   6.00000
REMARK    THIS IS A SIMULATION BOX
CRYST1   32.503   32.503   32.503  90.00  90.00  90.00 P 1           1
MODEL       16
ATOM      1  O1  THF     1      19.380  32.070  27.250  1.00  0.00
ATOM      2  CA1 THF     1      19.090  30.880  26.540  1.00  0.00
ATOM      3  CA2 THF     1      19.460  33.170  26.380  1.00  0.00
ATOM      4  CB1 THF     1      18.940  31.200  25.060  1.00  0.00
ATOM      5  CB2 THF     1      19.190  32.710  24.950  1.00  0.00
ATOM      6  T21 THF     1      18.950  30.900  25.460  1.00  0.00
ATOM      7  T22 THF     1      19.290  33.050  25.300  1.00  0.00
ATOM      8  O1  THF     2      31.020  24.150   7.910  1.00  0.00
ATOM      9  CA1 THF     2      31.310  24.450   9.250  1.00  0.00
ATOM     10  CA2 THF     2      31.450  22.850   7.570  1.00  0.00
ATOM     11  CB1 THF     2      32.010  23.250   9.900  1.00  0.00
ATOM     12  CB2 THF     2      32.100  22.200   8.790  1.00  0.00
ATOM     13  T21 THF     2      31.810  23.710   9.890  1.00  0.00
ATOM     14  T22 THF     2      31.950  22.220   8.310  1.00  0.00


so in this case the first 7 lines are data I don't need while the next 14 lines are important data I want to read. I know I can combine Skip and Readlist

str = OpenRead["file"];
Skip[str, String, 7]; coord = ReadList[str, {Word, Number, Word, Word}~Join~Table[Number, {6}], 14]


How can I do this repeatedly to read say 500 such data (i.e. ATOM 1 - 14, 500 times). I will like to write a function that can accept the number of such data I want to read.

This entry is related but not exactly what I want. If there's already a similar question I'll be more than happy to delete this but I couldn't find any that suits my purpose. Please note that these files are really huge ~1GB.

-
Since all those lines start with "ATOM", how about FindList? – Xerxes Mar 7 '13 at 0:51
@Xerxes, thanks for that idea, but some of the files don't all start with ATOM in those important rows. Some just have an index number. – RunnyKine Mar 7 '13 at 1:00
If you know some java you could probably edit Leonid's answer to my question that you linked to do what you want with some Buffered String readers. See: stackoverflow.com/questions/3806062/… – s0rce Mar 7 '13 at 2:00
@s0rce, you're probably right but I know absolutely zero Java. – RunnyKine Mar 7 '13 at 2:07
@RunnyKine I re-familiarize myself with it every morning by getting out my mug, and brewing a cup, or four. (Yes, the pun was intended.) – rcollyer Mar 7 '13 at 3:47

## 5 Answers

Simon's method is likely one of the fastest approaches. (Incidentally I believe it requires a little tweaking to work on your data, because there is extra space after the last number read that also needs to be skipped.)

Nevertheless I find it less that optimally convenient because you either need to know in advance the length of your file or check for EndOfFile within the loop. Another approach that is similarly fast is to describe an entire cycle within the specification for ReadList. This means we need a way to skip elements, and to do that I will use the method librik recommended here. (See this answer for more.)

I'll name my function skip:

skip[__String] = Sequence[];


I build the spec (the final skip[String] is the tweak I mentioned):

line = {Word, Number, Word, Word} ~Join~ Table[Number, {6}];
drop = skip @@ Table[Record, {7}];
spec = Join[{drop}, line ~Table~ {14}, {skip[String]}];


Which looks like:

Short[spec, 2]

{skip[Record,Record,Record,Record,Record,Record,Record],
{Word,Number,Word,Word,Number,Number,Number,Number,Number,Number},<<13>>,
skip[String]}


Now we pass this spec to ReadList which takes care of the rest:

ReadList["RKdata.txt", spec]


The third parameter of ReadList can be used to control how many blocks are read, e.g. ReadList["RKdata.txt", spec, 3] will read three blocks of 14 records.

You can use Join @@ or ~Flatten~ 1 if you want all of your records in a single list.

-
It took me longer than I care to admit to understand what is happening, and neither of the links explain "why" it works. It seems, I always assumed that the second argument of Read must be a list. The fact that it doesn't have to be any expression with Number, String, etc. used as placeholders makes it much more useful than I had previously understood it to be. +1 – rcollyer Mar 8 '13 at 19:31
@rcollyer The documentation for Read states: "You can use Read to get objects to insert into any expression structure, not necessarily a list. Example: Read[stream,Hold[Expression]] gets an expression and places it inside Hold." Nevertheless it wasn't until librik's comment (my first link) that I realized just what that implied. – Mr.Wizard Mar 8 '13 at 22:33
I'm aware of that, now. But, it is a single line buried in the "Details," so I've missed it up to this point. And, the implication isn't made clear until it is seen in action. Definitely something I have to play with. – rcollyer Mar 8 '13 at 23:43
This is brilliant. Thanks. – RunnyKine Mar 11 '13 at 0:57
@RunnyKine Happy to help. :-) – Mr.Wizard Mar 11 '13 at 1:02

Since you know how many lines to skip and how many to read, you can just wrap your code in a Table:

lineformat = {Word, Number, Word, Word}~Join~Table[Number, {6}];

str = OpenRead["file"];
coords = Table[Skip[str, String, 7]; ReadList[str, lineformat, 14], {500}];
Close[str];


This will of course create a nested list of dimensions {500,14,12}. If you need to remove the outer level of nesting, so that coords has dimensions {7000,12}, you can use either:

coords = Flatten[coords, 1];


or

coords = Join @@ coords;

-

This isn't a purely Mathematica-based answer, but an alternative approach might be to use something like the command-line utility sed to access the lines:

getLines[start_, end_, file_] :=
ReadList[
ToString@StringForm["!sed -n ', p' < ", start, end, file],
String,
RecordLists -> True]


which you can use like this:

getLines[5000, 5010, "/private/var/log/appfirewall.log"]


to produce 11 (:-)) lines of warnings from the log file.

-
Or awk for something a bit more complex. Even if it sounds heretic, this is probably the best approach for large files ... +1 – Szabolcs Mar 7 '13 at 17:39
@szabolcs oh yes, awk would be good too. For a 1GB text file, sed can retrieve lines in a second or so. – cormullion Mar 7 '13 at 18:45
ReadList["!perl -pn <perlcode>",...]


or

ReadList["!grep <regex>",...]


for example

ReadList["!grep '^[A-Z].*[0-9]\{1,2\}.*THF.*[0-9].*' textsample.txt",String]


In the perl or grep you put either pattern matcher / regex if you don't know how many lines, or you put something based on counts, if you DO know how many lines you want to take/skip. The former, regex, might be safer.

My grep example assumes that the first colum consists of upper-case letters, then come one or two integers, further right comes THF, and then come numbers. If this assumption is not correct, it will likely break.

Note that ReadList is, unlike Import, extremely fast and uses native i/o methods. (It has the same speed as OpenRead which is the fastest way to read data from a file into M).

This method also has the advantage that you only read your relevant extraction into the kernel. Stuff you don't want doesn't even make it into the kernel, further improving your efficiency. And perl is JIT-compiled and therefore runs with C-compiled speed. You could also use python if you prefer python over perl. grep is also machine-code speed.

-
Would you add specific PERL code for this particular data file? I'd like to vote for this method but I think that should be included. (BTW please note that not all lines start with "ATOM" as stated in the comments below the question.) – Mr.Wizard Mar 7 '13 at 23:57
@Mr. Wizard, that would depend on certain assumptions about the data that would be risky for me to make. That's why I intentionally left that unspecified. For example, if I assume ten columns, we'd also capture the rows that start with CRYST1. If I could assume the last five columns are reals and columns 1 and 4 are only uppercase characters, we'd come closer. – Andreas Lauschke Mar 8 '13 at 0:24
@Mr.Wizard, I've updated it, based on assumptions on the pattern matching with regex. As I said, it's a bit risky, because it's not further specified, we only have a short sample. – Andreas Lauschke Mar 8 '13 at 15:33
+1 for jumping through my hoops. :o) – Mr.Wizard Mar 8 '13 at 22:52

You want to isolate entries starting with "ATOM"?

I would read this as a Table and winnow entries using Select.

idata = Import["File", "Table"];
atom = Select[idata, SameQ[#[[1]], "ATOM"] &];

-
I would rather not use Import, since these files are really huge ~1GB. – RunnyKine Mar 7 '13 at 0:57
Alas, currently not an option for the files of ~ Gb size, see e.g. this discussion – Leonid Shifrin Mar 7 '13 at 1:07