The structural stability and conformational flexibility of two phycocyanobilin (PCB) binding bacteriaphytochromes, namely, Cph 1 from Synechocystis and the GAF domain of SyB from Synechococcus, were studied using all-atoms molecular dynamics simulations techniques. In order to involve the tetrapyrrole cofactor in the simulation, new empirical force field parameters were developed for PCB which are compatible with the CHARMM22 force field for proteins. Special emphasis was made in understanding the conflicting NMR-based structures recently obtained for the two parent states, Pr and Pfr, of SyB(GAF) regarding the highly distorted cofactor conformation which is in contrast to all crystallographic measurements on other phytochrome species. The 25 ns all-atoms MD simulation of Cph 1 in the Pr state and SyB(GAF) in its two parent states supports the picture of a relatively planar conformation of rings A, B, and C and ring D slightly twisted out of the ABC plane, in good agreement with crystallographic and spectroscopic experiments. SyB(GAF) converges to two very similar structures for the Pr state and also for the Pfr state with a ZZZssa conformation in disagreement with the ZZE configuration observed experimentally for other phytochromes. The failure in the prediction of the SyB(GAF) Pfr geometry is a consequence of deformed initial structure. Furthermore, in contrast to the results obtained for the SyB(GAF), the MD simulations showed a very stable Cph 1 structure where an important hydrogen bond and water network in the chromophore binding pocket could be identified. We did not find evidence for structural heterogeneity either for Cph 1 or SyB(GAF) at least on the nanosecond time scale.