Ignore:
Timestamp:
Jul 6, 2011 5:23:09 PM (10 years ago)
Author:
toby
Message:

a bit of cleanup

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/sandbox/readexp.tcl

    r1115 r1152  
    19131913            }
    19141914            abscor2-set {
    1915                 # can't use validreal as the decimal must be in col 20
    1916                 if {[catch {
    1917                     if {abs($value) < 99.99 && abs($value) > 1.e-4} {
    1918                         set tmp [format "%15.10f" $value]
    1919                         # make a final check of decimal
    1920                         if {[string range $tmp 4 4] != "."} {
    1921                             set tmp [format "%15.6E" $value]
    1922                         }
    1923                     } else {
    1924                         set tmp [format "%15.6E" $value]
     1915                # this must have a decimal as the 5th character, so that we end up with a
     1916                # decimal point in column 20.
     1917                set tmp $value
     1918                if ![validreal tmp 12 7] {return 0}
     1919                set pos [string first "." $tmp]
     1920                while {$pos < 4} {
     1921                    set tmp " $tmp"
     1922                    set pos [string first "." $tmp]
     1923                }
     1924                if {$pos == 4} {
     1925                    setexp "${key}ABSCOR" $tmp 16 15
     1926                    return
     1927                }
     1928                catch {
     1929                    set tmp [format "%12.6E" $value]
     1930                    set pos [string first "." $tmp]
     1931                    while {$pos < 4} {
     1932                        set tmp " $tmp"
     1933                        set pos [string first "." $tmp]
    19251934                    }
    1926                 }]} {return 0}
    1927                 setexp "${key}ABSCOR" $tmp 16 15
     1935                    if {$pos == 4} {
     1936                        setexp "${key}ABSCOR" $tmp 16 15
     1937                        return
     1938                    }
     1939                }
     1940                return 0
    19281941            }
    19291942            abstype-get {
     
    21692182
    21702183proc atom_constraint_read {phase} {
    2171 
    2172          set fix_list ""
    2173          foreach k {1 2 3 4 5 6 7 8 9} {
    2174                   set key [format "LEQV HOLD%1d%2d" $phase $k]
    2175                   set line [readexp $key]
    2176                   foreach j {2 10 18 26 34 42 50 58} {
    2177                            set fix_param [string range $line $j [expr $j+7]]
    2178                            if {[string trim $fix_param] == ""} {return $fix_list}
    2179                            lappend fix_list $fix_param
    2180                   }
    2181          }
     2184    set fix_list ""
     2185    foreach k {1 2 3 4 5 6 7 8 9} {
     2186        set key [format "LEQV HOLD%1d%2d" $phase $k]
     2187        set line [readexp $key]
     2188        foreach j {2 10 18 26 34 42 50 58} {
     2189            set fix_param [string range $line $j [expr $j+7]]
     2190            if {[string trim $fix_param] == ""} {return $fix_list}
     2191            lappend fix_list $fix_param
     2192        }
     2193    }
    21822194}
    21832195
    21842196# load all atom constraints into global array fix_param
    21852197proc atom_constraint_load { } {
    2186      catch {unset ::fix_param}
    2187      foreach i $::expmap(phaselist) {
    2188              set temp [atom_constraint_read $i]
    2189              foreach j $temp {
    2190                      set atomnum [string trim [string range $j 2 3]]
    2191                      set param [string trim [string range $j 4 6]]
    2192                      set ::fix_param($i,$atomnum,$param) 1
    2193 
    2194                      }
    2195      }
     2198    catch {unset ::fix_param}
     2199    foreach i $::expmap(phaselist) {
     2200        set temp [atom_constraint_read $i]
     2201        foreach j $temp {
     2202            set atomnum [string trim [string range $j 2 3]]
     2203            set param [string trim [string range $j 4 6]]
     2204            set ::fix_param($i,$atomnum,$param) 1   
     2205        }
     2206    }
    21962207}
    21972208
    21982209proc atom_constraint_write {phase fix_list} {
    2199 
    2200      foreach key [array names ::exparray "LEQV HOLD$phase*"] {
    2201              delexp $key
    2202      }
    2203      set k 0
    2204      set j 1
    2205      set line ""
    2206      foreach fix $fix_list {
    2207              incr k
    2208              append line $fix
    2209              if {$k == 8} {
    2210                 set key [format "LEQV HOLD%1d%2d" $phase $j]
    2211                 makeexprec $key
    2212                 setexp $key $line 3 [expr ($k * 8) + 2]
    2213                 set k 0
    2214                 incr j
    2215                 set line ""
    2216              }
    2217      }
    2218      if {$line != ""} {
     2210    foreach key [array names ::exparray "LEQV HOLD$phase*"] {
     2211        delexp $key
     2212    }
     2213    set k 0
     2214    set j 1
     2215    set line ""
     2216    foreach fix $fix_list {
     2217        incr k
     2218        append line $fix
     2219        if {$k == 8} {
     2220            set key [format "LEQV HOLD%1d%2d" $phase $j]
     2221            makeexprec $key
     2222            setexp $key $line 3 [expr ($k * 8) + 2]
     2223            set k 0
     2224            incr j
     2225            set line ""
     2226        }
     2227    }
     2228    if {$line != ""} {
    22192229        set key [format "LEQV HOLD%1d%2d" $phase $j]
    22202230        makeexprec $key
    22212231        setexp $key $line 3 [expr ($k * 8) + 2]
    2222      }
    2223 
    2224 }
    2225 
     2232    }   
     2233}
    22262234
    22272235
     
    35043512    }
    35053513}
     3514#============================================================================
     3515# rigid body EXP editing routines (to move into readexp.tcl)
     3516# RigidBodyList -- returns a list of the defined rigid body types
     3517# SetRigidBodyVar -- set variables and damping for rigid body type multipliers
     3518# ReadRigidBody  -- # of times a body is mapped, scaling factors, var #s & coordinates
     3519# RigidBodyMappingList - return a list instances where a RB is mapped in phase
     3520# RigidBodyEnableTLS -- Enable or Disable TLS use for a rigid body mapping
     3521# RigidBodySetTLS  -- change the TLS values for a rigid body mapping
     3522# RigidBodySetDamp -- change the damping values for a rigid body mapping
     3523# RigidBodyVary    -- set refinement variable numbers for a rigid body mapping
     3524# RigidBodyTLSVary -- set TLS refinement variable nums for a rigid body mapping
     3525# AddRigidBody -- defines a new rigid body type
     3526# DeleteRigidBody -- remove a rigid body definition
     3527# ReplaceRigidBody -- replaces a previous rigid body type
     3528# ReadRigidBodyMapping  -- get parameters for a rigid body mapping
     3529# MapRigidBody -- map a rigid body type into a phase
     3530# EditRigidBodyMapping -- change the parameters in a rigid body mapping
     3531# UnMapRigidBody --remove a rigid body constraint by removing a RB "instance"
     3532#----- note that these older routines should not be used ------
     3533# RigidBodyCount -- returns the number of defined rigid bodies (body types)
     3534#    use RigidBodyList instead
     3535# RigidBodyMappingCount -- # of times a rigid body is mapped in phase
     3536#    use RigidBodyMappingList instead
     3537#============================================================================
     3538# returns the number of defined rigid bodies
     3539proc RigidBodyCount {} {
     3540    set n [string trim [readexp "RGBD  NRBDS"]]
     3541    if {$n == ""} {
     3542        set n 0
     3543    }
     3544    return $n
     3545}
     3546
     3547# returns a list of the defined rigid body types
     3548proc RigidBodyList {} {
     3549    set n [string trim [readexp "RGBD  NRBDS"]]
     3550    if {$n == ""} {
     3551        set n 0
     3552    }
     3553    set rblist {}
     3554    foreach rbnum {1 2 3 4 5 6 7 8 9 10 11 12 13 14 15} {
     3555        set value $rbnum
     3556        validint value 2
     3557        set key "RGBD${value}"
     3558        if {[existsexp "$key NATR "]} {
     3559            lappend rblist $rbnum
     3560        }
     3561        if {[llength $rblist] == $n} break
     3562    }
     3563    return $rblist
     3564}
     3565
     3566# ReadRigidBody provides all information associated with a rigid body type
     3567#  rbnum is the rigid body type number
     3568# it returns two items:
     3569#   the number of times the rigid body is mapped
     3570#   a list containing an element for each scaling factor in rigid body #rbnum.
     3571# in each element there are four items:
     3572#    a multiplier value for the rigid body coordinates
     3573#    a damping value (0-9) for the refinement of the multiplier
     3574#    a variable number if the multiplier will be refined
     3575#    a list of cartesian coordinates coordinates
     3576# each cartesian coordinate contains 4 items: x,y,z and a label
     3577#  note that the label is present only when the RB is created in EXPGUI and is
     3578#  not used in GSAS.
     3579proc ReadRigidBody {rbnum} {
     3580    if {[lsearch [RigidBodyList] $rbnum] == -1} {
     3581        return ""
     3582    }
     3583    set value $rbnum
     3584    validint value 2
     3585    set key "RGBD${value}"
     3586    set n [string trim [string range [readexp "$key NATR"] 0 4]]
     3587    set used [string trim [string range [readexp "$key NBDS"] 0 4]]
     3588    set nmult [string trim [string range [readexp "$key NSMP"] 0 4]]
     3589    set out {}
     3590    for {set i 1} {$i <= $nmult} {incr i} {
     3591        set line [readexp "${key}${i}PARM"]
     3592        set mult [string trim [string range $line 0 9]]
     3593        set var [string trim [string range $line 10 14]]
     3594        set damp [string trim [string range $line 15 19]]
     3595        set coordlist {}
     3596        for {set j 1} {$j <= $n} {incr j} {
     3597            set value $j
     3598            validint value 3
     3599            set line [readexp "${key}${i}SC$value"]
     3600            set x [string trim [string range $line 0 9]]
     3601            set y [string trim [string range $line 10 19]]
     3602            set z [string trim [string range $line 20 29]]
     3603            set lbl [string trim [string range $line 30 39]]
     3604            lappend coordlist [list $x $y $z $lbl]
     3605        }
     3606        lappend out [list $mult $damp $var $coordlist]
     3607    }
     3608    return [list $used $out]
     3609}
     3610
     3611# SetRigidBodyVar
     3612#   rbnum is the rigid body type number
     3613#   varnumlist is a list of variable numbers
     3614#      note that if this list is shorter than the number of actual multipliers
     3615#      for the body, the unspecified variable will not be changed
     3616#   damplist   is a list of damping values (0-9)
     3617#      note that if the damplist is shorter than the number of actual multipliers
     3618#      the unspecified values are not changed
     3619#  SetRigidBodVar 2 {1 2 3} {}
     3620#       will vary the (first 3) translations in body #3 and will not change the
     3621#       damping values
     3622#  SetRigidBodVar 3 {} {0 0 0}
     3623#       will not change variable settings but will change the (first 3) damping values
     3624#  SetRigidBodVar 4 {11 11} {8 8}
     3625#      changes both variable numbers and damping at the same time
     3626# Nothing is returned
     3627proc SetRigidBodyVar {rbnum varnumlist damplist} {
     3628    if {[lsearch [RigidBodyList] $rbnum] == -1} {
     3629        return ""
     3630    }
     3631    set value $rbnum
     3632    validint value 2
     3633    set key "RGBD${value}"
     3634    set nmult [string trim [string range [readexp "$key NSMP"] 0 4]]
     3635    for {set i 1} {$i <= $nmult} {incr i} {
     3636        set j $i
     3637        incr j -1
     3638        set var [lindex $varnumlist $j]
     3639        if {$var != ""} {
     3640            validint var 5
     3641            setexp "${key}${i}PARM" $var 11 15
     3642        }
     3643        set damp [lindex $damplist $j]
     3644        if {$damp != ""} {
     3645            if {$damp > 9} {set damp 9}
     3646            if {$damp < 0} {set damp 0}
     3647            validint damp 5
     3648        }
     3649        setexp "${key}${i}PARM" $damp 16 20
     3650    }
     3651}
     3652
     3653
     3654# return the number of times rigid body $bodytyp is mapped in phase $phase
     3655proc RigidBodyMappingCount {phase bodytyp} {
     3656    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1]"
     3657    if {! [existsexp "$key  NBDS"]} {return 0}
     3658    set n [string trim [readexp "$key  NBDS"]]
     3659    if {$n == ""} {
     3660        set n 0
     3661    }
     3662    return $n
     3663}
     3664# return a list of the instances where rigid body $bodytyp is mapped in phase $phase
     3665proc RigidBodyMappingList {phase bodytyp} {
     3666    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1]"
     3667    if {! [existsexp "$key  NBDS"]} {return {}}
     3668    set n [string trim [readexp "$key  NBDS"]]
     3669    if {$n == ""} {
     3670        set n 0
     3671    }
     3672    set rblist {}
     3673    foreach rbnum {1 2 3 4 5 6 7 8 9 10 11 12 13 14 15} {
     3674        set value $rbnum
     3675        validint value 2
     3676        set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $rbnum 1]"
     3677        if {[existsexp "$key  NDA"]} {
     3678            lappend rblist $rbnum
     3679        }
     3680        if {[llength $rblist] == $n} break
     3681    }
     3682    return $rblist
     3683}
     3684
     3685
     3686
     3687# reads rigid body mapping parameters for phase ($phase), body type # ($bodytyp) and instance # ($num)
     3688# returns a list of items (most lists) as follows:
     3689#   1) sequence # of first atom in body
     3690#   2) origin of body in fractional coordinates (3 elements)
     3691#   3) Euler angles as 6 pairs of numbers (see below)
     3692#   4) variable numbers for the 9 position variables (origin followed by rotations)
     3693#   5) damping vals for the 9 position variables (origin followed by rotations)
     3694#   6) the TLS values, in order below (empty list if TLS is not in use)
     3695#   7) the variable numbers for each TLS values, in order below (or empty)
     3696#   8) three damping values for the T, L and S terms.
     3697# returns an empty list if no such body exists.
     3698#
     3699# Euler angles are a list of axes and angles to rotate:
     3700#   { {axis1 angle1} {axis2 angle2} ...}
     3701# where axis1,... can be 1, 2 or 3 corresponding to the cartesian X, Y or Z axes
     3702#
     3703# The 20 TLS terms are ordered:
     3704#    T11, T22, T33, T12, T13, T23
     3705#    L11, L22, L33, L12, L13, L23
     3706#    S12, S13, S21, S23, S31, S32, SAA, SBB
     3707#
     3708proc ReadRigidBodyMapping {phase bodytyp num} {
     3709    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
     3710        return ""
     3711    }
     3712    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
     3713    set first [string trim [string range [readexp "$key  NDA"] 0 4]]
     3714    set line [readexp "$key BDFL"]
     3715    set varlist {}
     3716    set damplist {}
     3717    foreach i {0 1 2 3 4 5 6 7 8} {
     3718        lappend varlist [string trim [string range $line [expr {5*$i}] [expr {4 + 5*$i}] ]]
     3719        lappend damplist [string trim [string range $line [expr {45 + $i}] [expr {45 + $i}] ]]
     3720    }
     3721    set TLSdamplist {}
     3722    foreach i {54 55 56} {
     3723        lappend TLSdamplist [string trim [string range $line $i $i ]]
     3724    }
     3725    set line [readexp "${key} BDLC"]
     3726    set x [string trim [string range $line 0 9]]
     3727    set y [string trim [string range $line 10 19]]
     3728    set z [string trim [string range $line 20 29]]
     3729    set origin [list $x $y $z]
     3730    set line [readexp "${key} BDOR"]
     3731    set rotations {}
     3732    foreach i {0 10 20 30 40 50} {
     3733        set angle [string trim [string range $line $i [expr {$i+7}]]]
     3734        set axis [string trim [string range $line [expr {$i+8}] [expr {$i+9}]]]
     3735        lappend rotations [list $angle $axis]
     3736    }
     3737    set TLS [string trim [string range [readexp "${key} LSTF"] 0 4]]
     3738    set tlsvars {}
     3739    set tlsvals {}
     3740    if {$TLS != 0} {
     3741        set line [readexp "${key}TLSF1"]
     3742        for {set j 0} {$j < 20} {incr j} {
     3743            set var [string trim [string range $line [expr {3*$j}] [expr {3*$j+2}]]]
     3744            if {$var == ""} {set var 0}
     3745            lappend tlsvars $var
     3746        }
     3747        for {set j 0} {$j < 20} {incr j} {
     3748            set i 0
     3749            if {$j == 0} {
     3750                set i 1
     3751            } elseif {$j == 8} {
     3752                set i 2
     3753            } elseif {$j == 16} {
     3754                set i 3
     3755            }
     3756            if {$i != 0} {
     3757                set line [readexp "${key}TLSP$i"]
     3758                set i 0
     3759                set j1 0
     3760                set j2 7
     3761            } else {
     3762                incr j1 8
     3763                incr j2 8
     3764            }
     3765            set val [string trim [string range $line $j1 $j2]]
     3766            if {$val == ""} {set val 0}
     3767            lappend tlsvals $val
     3768        }
     3769    }
     3770    return [list $first $origin $rotations $varlist $damplist $tlsvals $tlsvars $TLSdamplist]
     3771}
     3772
     3773# Control TLS representation for phase, body # and instance number of a Rigid body mapping
     3774#   for mapping with phase ($phase), body type # ($bodytyp) and instance # ($num)
     3775# Enable TLS use if TLS is non-zero (true). Disable if zero
     3776proc RigidBodyEnableTLS {phase bodytyp num TLS} {
     3777    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
     3778        return ""
     3779    }
     3780    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
     3781    if {$TLS} {
     3782        setexp "${key} LSTF" [format "%5d" 1] 1 5
     3783        if {![existsexp "${key}TLSF1"]} {makeexprec "${key}TLSF1"}
     3784        if {![existsexp "${key}TLSP1"]} {
     3785            makeexprec "${key}TLSP1"
     3786            set str {}
     3787            foreach v {.01 .01 .01 0 0 0 0 0} d {4 4 4 4 4 4 2 2} {
     3788                validreal v 8 $d
     3789                append str $v
     3790            }
     3791            setexp "${key}TLSP1" $str 1 64
     3792        }
     3793        if {![existsexp "${key}TLSP2"]} {
     3794            makeexprec "${key}TLSP2"
     3795            set str {}
     3796            set v 0
     3797            foreach d {2 2 2 2 4 4 4 4} {
     3798                validreal v 8 $d
     3799                append str $v
     3800            }
     3801            setexp "${key}TLSP2" $str 1 64
     3802        }
     3803        if {![existsexp "${key}TLSP3"]} {
     3804            makeexprec "${key}TLSP3"
     3805            set str {}
     3806            set v 0
     3807            foreach d {4 4 4 4} {
     3808                validreal v 8 $d
     3809                append str $v
     3810            }
     3811            setexp "${key}TLSP3" $str 1 64
     3812        }
     3813    } else {
     3814        setexp "${key} LSTF" [format "%5d" 0] 1 5
     3815    }
     3816    return 1
     3817}
     3818
     3819# Control the TLS values for Rigid body mapping for mapping with
     3820#    phase ($phase), body type # ($bodytyp) and instance # ($num)
     3821# set the 20 TLS values to the values in TLSvals
     3822# There must be exactly 20 TLS terms, which are ordered:
     3823#    T11, T22, T33, T12, T13, T23
     3824#    L11, L22, L33, L12, L13, L23
     3825#    S12, S13, S21, S23, S31, S32, SAA, SBB
     3826proc RigidBodySetTLS {phase bodytyp num TLSvals} {
     3827    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
     3828        return ""
     3829    }
     3830    if {[llength $TLSvals] != 20} {return ""}
     3831    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
     3832    set TLS [string trim [string range [readexp "${key} LSTF"] 0 4]]
     3833    if {$TLS == 0} {return ""}
     3834    if {![existsexp "${key}TLSF1"]} {makeexprec "${key}TLSF1"}
     3835    foreach n {1 2 3} {
     3836        if {![existsexp "${key}TLSP$n"]} {makeexprec "${key}TLSP$n"}
     3837    }
     3838    set str {}
     3839    set n 1
     3840    set i 0
     3841    foreach v $TLSvals d {4 4 4 4 4 4 2 2 2 2 2 2 4 4 4 4 4 4 4 4} {
     3842        incr i
     3843        validreal v 8 $d
     3844        append str $v
     3845        if {$i == 8} {
     3846            set i 0
     3847            setexp "${key}TLSP$n" $str 1 64
     3848            incr n
     3849            set str {}
     3850        }
     3851    }
     3852    setexp "${key}TLSP$n" $str 1 64
     3853    return 1
     3854}
     3855
     3856# set damping values for a Rigid body mapping
     3857#   for mapping with phase ($phase), body type # ($bodytyp) and instance # ($num)
     3858# there must be 9 damping values in RBdamp for the 9 position variables (origin followed by rotations)
     3859# Use of TLSdamp is optional, but to be used, TLS representation must be enabled and there must be
     3860# three damping terms (for all T terms; for all L terms and for all S terms)
     3861proc RigidBodySetDamp {phase bodytyp num RBdamp "TLSdamp {}"} {
     3862    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
     3863        return ""
     3864    }
     3865    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
     3866    if {[llength $RBdamp] != 9} {return ""}
     3867    set str {}
     3868    foreach v $RBdamp {
     3869        if {[validint v 1] != 1} {set v " "}
     3870        append str $v
     3871    }
     3872    setexp "$key BDFL" $str 46 9
     3873    set TLS [string trim [string range [readexp "${key} LSTF"] 0 4]]
     3874    if {$TLS != 0 &&  [llength $TLSdamp] == 3} {
     3875        set str {}
     3876        foreach v $TLSdamp {
     3877        if {[validint v 1] != 1} {set v " "}
     3878            append str $v
     3879        }
     3880        setexp "$key BDFL" $str 55 3
     3881    }
     3882    return 1
     3883}
     3884
     3885# set refinement variable numbers for a Rigid body mapping
     3886#   for mapping with phase ($phase), body type # ($bodytyp) and instance # ($num)
     3887# there must be 9 variable values in RBvar for the 9 position variables (origin followed by rotations)
     3888# note that the variable values should be unique integers
     3889proc RigidBodyVary {phase bodytyp num RBvar} {
     3890    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
     3891        return ""
     3892    }
     3893    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
     3894    if {[llength $RBvar] != 9} {return ""}
     3895    set str {}
     3896    foreach v $RBvar {
     3897        if {[validint v 5] != 1} {set v " "}
     3898        append str $v
     3899    }
     3900    setexp "$key BDFL" $str 1 45   
     3901}
     3902
     3903# set TLS refinement variable numbers for a Rigid body mapping
     3904#   for mapping with phase ($phase), body type # ($bodytyp) and instance # ($num)
     3905# there must be 20 variable values in TLSvar for the 20 parameters:
     3906#    T11, T22, T33, T12, T13, T23
     3907#    L11, L22, L33, L12, L13, L23
     3908#    S12, S13, S21, S23, S31, S32, SAA, SBB
     3909# note that the variable values should be unique integers
     3910proc RigidBodyTLSVary {phase bodytyp num TLSvar} {
     3911    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
     3912        return ""
     3913    }
     3914    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
     3915    if {[llength $TLSvar] != 20} {return ""}
     3916    set TLS [string trim [string range [readexp "${key} LSTF"] 0 4]]
     3917    if {$TLS == 0} {return ""}
     3918    set str {}
     3919    foreach v $TLSvar {
     3920        if {[validint v 3] != 1} {set v " "}
     3921        append str $v
     3922    }
     3923    setexp "${key}TLSF1" $str 1 60
     3924
     3925# AddRigidBody: add a new rigid body definition into the .EXP file
     3926# arguments are:
     3927#   multlist: defines a list of multipliers for each set of coordinates. In the
     3928#             simplest case this will be {1}
     3929#   coordlist: a nested list of coordinates such as { { {0 0 0} {.1 .1 .1} {.2 .2 .2} } }
     3930# note that when the length of multlist > 1 then coordlist must have the same length.
     3931# for input where
     3932#     multlist = {s1 s2} and
     3933#     coordlist = { { {0 0 0} {1 1 0} {.0 .0 .0} ...}
     3934#                     {0 0 0} {1 1 0} {2 1 2} ...}
     3935#                 }
     3936# the cartesian coordinates are defined from the input as
     3937#    atom 1 = s1 * (0,0,0) + s2*(0,0,0) [ = (0,0,0)]
     3938#    atom 2 = s1 * (1,1,0) + s2*(1,1,0) [ = (s1+s2) * (1,1,0)]
     3939#    atom 3 = s1 * (0,0,0) + s2*(2,1,2) [ = s2 * (2,1,2)]
     3940#    ...
     3941# Returns the number of the rigid body that has been created
     3942proc AddRigidBody {multlist coordlist} {
     3943    # find the first unused body #
     3944    foreach rbnum {1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16} {
     3945        set value $rbnum
     3946        validint value 2
     3947        set key "RGBD${value}"
     3948        if {! [existsexp "$key NATR "]} {break}
     3949    }
     3950    # did we go too far?
     3951    if {$rbnum == 16} {return ""}
     3952    # increment the RB counter
     3953    set n [string trim [readexp "RGBD  NRBDS"]]
     3954    if {$n == ""} {
     3955        makeexprec "RGBD  NRBDS"
     3956        set n 0
     3957    }
     3958    incr n
     3959    validint n 5
     3960    setexp "RGBD  NRBDS" $n 1 5
     3961    SetRigidBody $rbnum $multlist $coordlist
     3962    return $rbnum
     3963}
     3964
     3965# DeleteRigidBody: remove a rigid body definition from the .EXP file
     3966# The body may not be mapped. I am not sure if GSAS allows more than 9 bodies,
     3967# but if it does, the simplifed approach used here will fail, so this
     3968# is not allowed.
     3969# Input:
     3970#   Rigid body number
     3971# Returns:
     3972#   1 on success
     3973#   -1 if the body number is 11 or greater
     3974#   -2 if the body is mapped
     3975#   -3 if the body is not defined
     3976proc DeleteRigidBody {rbnum} {
     3977    # can't delete bodies with numbers higher than 10, since the key prefix
     3978    # (RGBD11... will overlap with rigid body instance records, which would be
     3979    # deleted below
     3980    if {$rbnum > 10} {
     3981        return -1
     3982    }
     3983    set value $rbnum
     3984    validint value 2
     3985    set key "RGBD${value}"
     3986    if {![existsexp "$key NATR "]} {
     3987        return -2
     3988    }
     3989    # make sure the body is not mapped
     3990    if {[string trim [string range [readexp "$key NBDS"] 0 4]] != 0} {
     3991        return -3
     3992    }
     3993    # delete the records starting with "RGBD x" or "RGBD10"
     3994    foreach key [array names ::exparray "${key}*"] {
     3995        #puts $key
     3996        delexp $key
     3997    }
     3998    # decrement the RB counter
     3999    set n [string trim [readexp "RGBD  NRBDS"]]
     4000    if {$n == ""} {
     4001        set n 0
     4002    }
     4003    incr n -1
     4004    validint n 5
     4005    if {$n > 0} {
     4006        setexp "RGBD  NRBDS" $n 1 5
     4007    } else {
     4008        delexp "RGBD  NRBDS"
     4009    }
     4010    return 1
     4011}
     4012
     4013# ReplaceRigidBody: replace all the information for rigid body #rbnum
     4014# Works the sames as AddRigidBody (see above) except that the rigid body is replaced rather
     4015# than added.
     4016# Note that count of the # of times the body is used is preserved
     4017proc ReplaceRigidBody {rbnum multlist coordlist {varlist ""} {damplist ""}} {
     4018    set value $rbnum
     4019    validint value 2
     4020    set key "RGBD${value}"
     4021    set line [readexp "$key NBDS"]
     4022    foreach key [array names ::exparray "${key}*"] {
     4023        #puts $key
     4024        delexp $key
     4025    }
     4026    SetRigidBody $rbnum $multlist $coordlist $varlist $damplist
     4027    setexp "$key NBDS" $line 1 68
     4028}
     4029
     4030# Edit the parameters for rigid body #rbnum
     4031# (normally called from ReplaceRigidBody or AddRigidBody)
     4032proc SetRigidBody {rbnum multlist coordlist {varlist ""} {damplist ""}} {
     4033    set value $rbnum
     4034    validint value 2
     4035    set key "RGBD${value}"
     4036    # number of atoms
     4037    set value [llength [lindex $coordlist 0]]
     4038    validint value 5
     4039    makeexprec "$key NATR"
     4040    setexp "$key NATR" $value 1 5
     4041    # number of times used
     4042    set value 0
     4043    validint value 5
     4044    makeexprec "$key NBDS"
     4045    setexp "$key NBDS" $value 1 5
     4046    # number of coordinate matrices
     4047    set value [llength $multlist]
     4048    validint value 5
     4049    makeexprec "$key NSMP"
     4050    setexp "$key NSMP" $value 1 5
     4051    set i 0
     4052    foreach mult $multlist coords $coordlist {
     4053        set var [lindex $varlist $i]
     4054        if {$var == ""} {set var 0}
     4055        set damp [lindex $damplist $i]
     4056        if {$damp == ""} {set damp 0}
     4057        incr i
     4058        makeexprec "${key}${i}PARM"
     4059        setexp "${key}${i}PARM" [format "%10.5f%5d%5d" $mult $var $damp] 1 20
     4060        set j 0
     4061        foreach item $coords {
     4062            #puts $item
     4063            incr j
     4064            set value $j
     4065            validint value 3
     4066            makeexprec "${key}${i}SC$value"
     4067            if {[llength $item] == 4} {
     4068                setexp "${key}${i}SC$value" [eval format "%10.6f%10.6f%10.6f%10s" $item] 1 40
     4069            } elseif {[llength $item] == 3} {
     4070                setexp "${key}${i}SC$value" [eval format "%10.6f%10.6f%10.6f" $item] 1 30
     4071            } else {
     4072                return -code 3 "Invalid number of coordinates"
     4073            }
     4074        }
     4075    }
     4076}
     4077
     4078# convert a decimal to the GSAS hex encoding with a field $digits long.
     4079proc ToHex {num digits} {
     4080    return [string toupper [format "%${digits}x" $num]]
     4081}
     4082
     4083# convert a GSAS hex encoding to a decimal integer
     4084proc FromHex {hex} {
     4085    return [scan $hex "%x"]
     4086}
     4087
     4088# MapRigidBody: define an "instance" of a rigid body: meaning that the coordinates
     4089# (and optionally U values) for a set of atoms will be generated from the rigid body
     4090# arguments:
     4091#   phase: phase number (1-9)
     4092#   bodytyp: number of rigid body (1-15) as returned from AddRigidBody
     4093#   firstatom: sequence number of the first atom in phase (note that atoms may
     4094#              not be numbered sequentially)
     4095#   position: list of three fractional coordinates for the origin of the rigid body coordinates
     4096#   angles: list of 3 angles to rotate the rigid body coordinates around x, y, z of the
     4097#           cartesian system before the body is translated to position.
     4098# returns the instance # (number of times body $bodytyp has been used in phase $phase)
     4099proc MapRigidBody {phase bodytyp firstatom position angles} {
     4100    # find the first unused body # for this phase & type
     4101    foreach rbnum {1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16} {
     4102        set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $rbnum 1]"
     4103        if {! [existsexp "$key  NDA"]} {break}
     4104    }
     4105    # did we go too far?
     4106    if {$rbnum == 16} {return ""}
     4107    # increment number of mapped bodies of this type overall
     4108    set value $bodytyp
     4109    validint value 2
     4110    set key "RGBD${value}"
     4111    set used [string trim [string range [readexp "$key NBDS"] 0 4]]
     4112    incr used
     4113    set value $used
     4114    validint value 5
     4115    setexp "$key NBDS" $value 1 5
     4116    # increment number of mapped bodies of this type in this phase
     4117    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1]"
     4118    if {[existsexp "$key  NBDS"]} {
     4119        set used [string trim [string range [readexp "$key  NBDS"] 0 4]]
     4120    } else {
     4121        makeexprec "$key  NBDS"
     4122        set used 0
     4123    }
     4124    incr used
     4125    set value $used
     4126    validint value 5
     4127    setexp "$key  NBDS" $value 1 5
     4128    # now write the mapping parameters
     4129    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $rbnum 1]"
     4130    set value $firstatom
     4131    validint value 5
     4132    makeexprec "$key  NDA"
     4133    setexp "$key  NDA" $value 1 5
     4134    set l1 {}
     4135    set l2 {}
     4136    for {set i 0} {$i < 9} {incr i} {
     4137        append l1 [format %5d 0]
     4138        append l2 [format %1d 0]
     4139    }
     4140    makeexprec "$key BDFL"
     4141    setexp "$key BDFL" $l1$l2 1 54
     4142    makeexprec "${key} BDLC"
     4143    setexp "${key} BDLC" [eval format "%10.6f%10.6f%10.6f" $position] 1 30
     4144    makeexprec "${key} BDOR"
     4145    set l1 {}
     4146    foreach val "$angles 0 0 0" dir "1 2 3 1 1 1" {
     4147        append l1 [format "%8.2f%2d" $val $dir]
     4148    }
     4149    setexp "${key} BDOR" $l1 1 60
     4150    makeexprec "${key} LSTF"
     4151    setexp "${key} LSTF" [format "%5d" 0] 1 5
     4152    return $rbnum
     4153}
     4154
     4155# EditRigidBodyMapping: edit parameters that define an "instance" of a rigid body (see MapRigidBody)
     4156# arguments:
     4157#   phase: phase number (1-9)
     4158#   bodytyp: number of rigid body (1-15) as returned from AddRigidBody
     4159#   bodynum: instance number, as returned by MapRigidBody
     4160#   position: list of three fractional coordinates for the origin of the rigid body coordinates
     4161#   angles: list of 3 angles to rotate the rigid body coordinates around x, y, z of the
     4162#           cartesian system before the body is translated to position.
     4163#
     4164proc EditRigidBodyMapping {phase bodytyp bodynum position angles} {
     4165    # number of bodies of this type in this phase
     4166    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $bodynum 1]"
     4167    setexp "${key} BDLC" [eval format "%10.6f%10.6f%10.6f" $position] 1 30
     4168    set l1 {}
     4169    foreach val "$angles 0 0 0" dir "1 2 3 1 1 1" {
     4170        append l1 [format "%8.2f%2d" $val $dir]
     4171    }
     4172    setexp "${key} BDOR" $l1 1 60
     4173}
     4174
     4175# UnMapRigidBody: remove a rigid body constraint by removing a RB "instance"
     4176# (undoes MapRigidBody)
     4177# arguments:
     4178#   phase: phase number (1-9)
     4179#   bodytyp: number of rigid body (1-15) as returned from AddRigidBody
     4180#   bodynum: instance number, as returned by MapRigidBody
     4181proc UnMapRigidBody {phase bodytyp mapnum} {
     4182    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $mapnum] == -1} {
     4183        return ""
     4184    }
     4185    # decrement number of mapped bodies of this type overall
     4186    set value $bodytyp
     4187    validint value 2
     4188    set key "RGBD${value}"
     4189    set used [string trim [string range [readexp "$key NBDS"] 0 4]]
     4190    incr used -1
     4191    set value $used
     4192    validint value 5
     4193    setexp "$key NBDS" $value 1 5
     4194    # decrement number of mapped bodies of this type in this phase
     4195    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1]"
     4196    if {[existsexp "$key  NBDS"]} {
     4197        set used [string trim [string range [readexp "$key  NBDS"] 0 4]]
     4198    } else {
     4199        set used 0
     4200    }
     4201    incr used -1
     4202    set value $used
     4203    validint value 5
     4204    if {$used > 0} {
     4205        setexp "$key  NBDS" $value 1 5
     4206    } else {
     4207        delexp "$key  NBDS"
     4208    }
     4209    # now delete the mapping parameter records
     4210    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $mapnum 1]"
     4211    foreach key [array names ::exparray "${key}*"] {
     4212        delexp $key
     4213    }
     4214    return $used
     4215}
     4216
Note: See TracChangeset for help on using the changeset viewer.