generates models based on number of connections
[openmx:openmx.git] / mxModelGen.R
1 # generates all models for given number of connections
2 # and writes them out as mxmodel objects
3
4 library(OpenMx)
5 dyn.load("/home/ac/skenny/perm.so")
6
7 ## should be like R_SWIFT_ARGS=3 3
8
9 ## get number of connections and array size (will be converted to nxn)
10
11 allinputs <- Sys.getenv("R_SWIFT_ARGS")
12
13 print(allinputs)
14
15 connections <- noquote(strsplit(allinputs," ")[[1]][1])
16
17 print(connections)
18
19 numcol        <- noquote(strsplit(allinputs," ")[[1]][2])
20
21 print(numcol)
22
23 rowcol = as.integer(numcol)
24 size = rowcol*rowcol
25
26
27 # creates the .adat files 
28
29 .C("run_perm", connections=as.integer(connections), size=as.integer(size))
30
31 # read files and generate associated ap.adat, s.dat and sp.dat files
32 # will do this in the c code
33
34 # read in the data files to generate models and write as objects
35
36 # dot is wildcard here, will need to adjust to it gets all 3 dat's
37 count = 0
38 for (i in list.files(pattern=".adat")){
39         model = MxModel()
40 #       model$A <- FullMatrix(rowcol, rowcol, free=TRUE)
41 #       model$S <- DiagMatrix(rowcol, rowcol) 
42         model$F <- IdenMatrix(rowcol, rowcol) # modifiable in future, for now no filter
43         model$I <- IdenMatrix(rowcol, rowcol)
44         model$cov <- MxAlgebra(model$F %&% (solve(model$I - model$A) %&% model$S))
45
46 # hard-coding for now
47 a = double(length = size)
48 apar = double(length = size)
49 xindex = array(c(1:rowcol))
50 yindex = array(c(1:rowcol))
51         apar = matrix(c(read.table(i)), nrow = rowcol, ncol = rowcol)
52         cnt = rowcol+1.0
53 for(x in xindex)
54 {
55  for(y in yindex)
56    {
57         if (apar[x,y] == 1){
58         apar[x,y] = cnt
59         cnt = cnt + 1}
60
61    }
62 }
63 print(apar)
64 aval = matrix(c(read.table(i)), nrow = rowcol, ncol = rowcol)
65
66 for(x in xindex)
67 {
68  for(y in yindex)
69    {
70         if (aval[x,y] == 1)
71         aval[x,y] = 0.75
72
73    }
74 }
75 print(aval)
76
77         
78         s = matrix(c(read.table("matrices/s.dat")), nrow = rowcol, ncol = rowcol)
79         spar = matrix(c(read.table("matrices/sp.dat")), nrow = rowcol, ncol = rowcol)
80         
81         model$A$values <- aval
82         model$A$parameters <- apar
83         model$S$values <- s
84         model$S$parameters <- spar
85         count = count+1
86 print(model$S)
87         save.Object(model,file=paste(count,"_model.rdata", sep=""))
88 }
89 system("tar -cf models.tar *rdata")